Derivative of a Wilson Line It's my first post on this website so please excuse any breaches of protocol that I'm unaware of. 
I've come across a formula for taking derivatives of Wilson Lines with respect to points on the path, but I'm having trouble deriving it. I found it in Equation 12 of https://arxiv.org/abs/1210.2581 , where it is presented without proof or citation. I therefore expect that it is an elementary proof that is very familiar to people. My attempt to prove it goes as follows:
$$W(x,y) \equiv \hat{P} \exp[ig\int A_\mu(s) d s^\mu] $$
Defines the path ordered Wilson line with starting point $y^\mu$ and endpoint $x^\mu$. Given a particular field $A_\mu$ I think of this as a functional of the path. Instead, if we can parametrize the path via $s^\mu(t)$, then the Wilson line can instead be though of a function of the endpoints:
$$W(x,y) = \hat{P} \exp[ig\int_{t_y}^{t_x} A_\mu(s(t))\frac{ ds^\mu}{dt} dt]$$
Where we've implied $s^\mu(t_y) = y^\mu$ and $s^\mu(t_x) = x^\mu$. I would like to get $\frac{\partial W(x,y)}{\partial x^\mu}$, but I only have $t_x$ as an accessible variable in the expression. So then I take $\frac{ \partial W(x,y)}{\partial t_x}$ instead:
$$\frac{ \partial W(x,y)}{\partial t_x} = igA_\mu(x) \frac{ ds^\mu}{d t}|_{t=t_x}W(x,y) $$
Then I think to myself that since $x^\mu$ is a function of $t_x$, maybe I can invert the function and use the chain rule:
$$\frac{ d}{dx^\mu} \overset?= \frac{dt_x}{dx^\mu}\frac{d}{dt_x}$$
(By this logic, it should really be $\frac{\partial t_x}{\partial x^\mu}$, but I'm thinking we can write $\frac{ d t_x}{d x^\mu} = (\frac{ d x^\mu}{d t_x})^{-1}$). Then, this $\frac{dt_x}{d x^\mu}$ looks like it will cancel the $\frac{ ds^\mu}{dt}|_{t=t_x}$ in the above expression since $s^\mu(t_x) = x^\mu$. So I write 
$$\frac{ \partial W(x,y)}{\partial x^\mu} = igA_\mu(x) W(x,y).$$
Note that if I took the derivative with respect to the start point $y^\mu$ instead, I would get a similar result but with an extra minus sign from the fundamental theorem of calculus:
$$\frac{ \partial W(x,y)}{\partial y^\mu} = -ig W(x,y) A_\mu(y).$$
(Throughout I have assumed that the chain rule maintains the ordering)
Then, if I want to extend this to taking derivatives with respect to any point on the path, not just an endpoint, I can use the fact that:
$$W(x,y) = W(x,z)W(z,y)$$
Straightforward application of the product rule and our above expressions then gives:
$$\frac{ \partial W(x,y)}{\partial z^\mu} = 0.$$
Which seems obvious because we've only allowed $z^\mu$ to move along the pre-specified path for the Wilson line. So the Wilson line itself does not change. The only way to get a nonzero derivative is to use the start point or the endpoint, because changing them changes the length of the entire path. However, in my understanding, this is actually incorrect. The result is not zero! Equation 12 of https://arxiv.org/abs/1210.2581 gives the result:
$$\frac{\partial W(x,y)}{\partial z^\mu} =  igW(x,s)A_\alpha(s) \frac{ \partial s^\alpha}{\partial z^\mu} W(s,y)|_{s=y}^{s=x}  + ig \int_{y}^{x} W(x,s) F_{\alpha \beta} (s) W(s,y) \frac{ \partial s^\alpha}{\partial z^\mu} ds^\beta. $$
Does anyone know how to derive this result? I stumbled on another paper: https://arxiv.org/abs/0805.3849 which seems to do the derivation, but only for an endpoint. They say that my result for the endpoint is the result of using something called the "Mandelstam condition." I'm not yet sure what this condition is, but hope to understand it soon! This paper also does a similar derivation to mine, and I can follow every step except one. The step in going from Eq. 21 to Eq. 22 is a complete mystery to me. If anyone can explain how that step was done, I would be extremely grateful! Thanks a bunch :)
 A: OP Here, I've figured out the issue. It really is a functional derivative. Below I include an excerpt from my own notes, in case another lost soul stumbles upon this page in the future. (Please excuse the changes in notation from my original post, I hope it won't be difficult to switch over, and please excuse the fact that some of this might be in excruciating detail - also certainly not the most efficient way of deriving the result)
$\newcommand{\D}{\partial}$
$\newcommand{\dd}{d}$
We've already encountered Wilson lines through a Feynman diagram analysis. It will turn out that they are a common object that arises when computing the propagator in the presence of a background field. In the most general form, Wilson lines are:
\begin{align}
V(z,y) \equiv \hat{P} \exp \left(ig \int_{y}^{z} A_\mu(s)\dd s^\mu\right).
\end{align}
$\hat{P}$ is the path ordering operator. The integral here is a line integral in spacetime from some starting point $y^\mu$ to some ending point $z^\mu$, along some path that we can parametrize via $s^\mu(t)$. 
\begin{align}
V(z,y) = \hat{P} \exp \left( ig \int_0^1 A_\mu(s(t)) \frac{ \dd s^\mu}{\dd t} \dd t\right).
\end{align}
We've defined $s^\mu(0) = y^\mu$ and $s^\mu(1) = z^\mu$. In order to see how this might appear as a solution to QCD equations of motion, we'll need to know how to take derivatives of Wilson lines. We assume here that the background field $A_\mu$ is a fixed function, and therefore the Wilson line is a functional of the path $s^\mu(t)$. We want to take a derivative with respect to spacetime coordinates $\frac{ \D }{\D x^\mu}$. Since $V(z,y)$ is a functional, we can think of it as a function of an infinite number of variables $s_i^\mu$. The label $i$ is purely symbolic, since it's really a function of $s(t)$ for every value of $t$. Using this symbolic notation, one can then use the chain rule:
\begin{align}
\frac{ \D V(z,y)}{\D x^\mu} = \frac{ \D V(z,y)}{\D s_i^\nu}\frac{ \D s_i^\nu}{\D x^\mu}.
\end{align}
Moving back to the functional notation, we can write $s_i^\nu = s^\nu(t)$ and then the implied sum over $i$ becomes an integral over $t$.
\begin{align} \label{Dform}
\frac{ \D V(z,y)}{\D x^\mu} = \int_0^1 \dd t \frac{ \delta V(z,y)}{\delta s^\nu(t)} \frac{ \D s^\nu}{\D x^\mu}.
\end{align}
Now it is clear that we cannot go further without figuring out the functional derivative of the Wilson line with respect to the path. Let's recall the method of finding the functional derivative. First, let's write $V(z,y)$ as a functional of $s^\nu(t)$.
\begin{align}
V[s(t)] = \hat{P} \exp\left( ig\int_0^1 \dd t A_\mu(s(t))\frac{ \dd s^\mu}{\dd t}\right).
\end{align}
Then let's make the transformation $s^\mu(t) \rightarrow s^\nu(t) + \delta s^\mu(t)$. We'll define $\delta V = V(s+\delta s) - V(s)$. Then, all we have to do is find $\delta V$ by simplifiying and ommitting any factors of $\mathcal{O}(\delta^2)$ or higher, and write the result as:
\begin{align} \label{funcderiv}
\delta V = \int_0^1 \dd t \frac{ \delta V}{\delta s^\mu(t)} \delta s^\mu(t),
\end{align}
and this will suffice to identify $ \frac{ \delta V}{\delta s^\mu(t)}$. So, let's begin as follows:
\begin{align}
\delta V = \hat{P}\exp\left(ig \int_0^1 \dd t A_\mu(s+\delta s)\left( \frac{ \dd s^\mu}{\dd t} + \frac{ \dd \delta s^\mu}{\dd t}\right)\right) - \hat{P} \exp\left(ig\int_0^1\dd t A_\mu(s) \frac{ \dd s^\mu}{\dd t}\right).
\end{align}
Taylor expanding $A_\mu$ in the first term to order $\delta$ gives
\begin{align}
\delta V = \hat{P}\exp \left( ig\int_0^1 \dd t \left( A_\mu(s) + \frac{ \D A_\mu}{\D s^\nu} \delta s^\nu\right)\left(\frac{ \dd s^\mu}{\dd t} + \frac{ \dd \delta s^\mu}{\dd t}\right)\right) -  \hat{P} \exp\left(ig\int_0^1\dd t A_\mu(s) \frac{ \dd s^\mu}{\dd t}\right).
\end{align}
Multiplying out the parentheses and ignoring the term of order $\delta^2$, we have
\begin{align}
\delta V = \hat{P}\exp \left( ig\int_0^1 \dd t \left( A_\mu(s) \frac{ \dd s^\mu}{\dd t} + A_\mu(s) \frac{ \dd \delta s^\mu}{\dd t} + \frac{ \D A_\mu}{\D s^\nu} \delta s^\nu \frac{ \dd s^\mu}{\dd t} \right)\right) -  \hat{P} \exp\left(ig\int_0^1\dd t A_\mu(s) \frac{ \dd s^\mu}{\dd t}\right).
\end{align}
Notice we can rewrite the $A_\mu(s) \frac{ \dd s^\mu}{\dd t}$ as $\frac{ \dd}{\dd t} \left( A_\mu(s) \delta s^\mu\right) - \frac{ \D A_\mu}{\D s^\nu} \frac{ \dd s^\nu}{\dd t} \delta s^\mu.$ [We used the chain rule to write $\frac{ \dd A_\mu}{\dd t} = \frac{ \D A_\mu}{\D s^\nu} \frac{ \dd s^\nu}{\dd t}$.] This is just integration by parts but keeping the boundary term. 
\begin{align} 
\delta V = \hat{P} \exp\left(ig \int_0^1 \dd t\left( A_\mu(s) \frac{ \dd s^\mu}{\dd t} + \frac{ \dd}{\dd t} (A_\mu(s) \delta s^\mu) + \left(\frac{ \D A_\mu}{\D s^\nu} - \frac{ \D A_\nu}{\D s^\mu}\right) \delta s^\nu \frac{ \dd s^\mu}{\dd t}\right)\right) - \hat{P} \exp\left(ig\int_0^1\dd t A_\mu(s) \frac{ \dd s^\mu}{\dd t}\right).
\end{align}
Now, notice that the path ordering operator will enfore ordering. Thus all matrices inside the exponent will commute. So we can insert commutators without changing the result, since any commutator is zero. Thus let's add an $ig[A_\mu,A_\nu]$ to complete the field strength tensor $F_{\nu\mu}$.
\begin{align} 
\delta V = \hat{P} \exp\left(ig \int_0^1 \dd t\left( A_\mu(s) \frac{ \dd s^\mu}{\dd t} + \frac{ \dd}{\dd t} (A_\mu(s) \delta s^\mu) + F_{\mu\nu} \delta s^\mu \frac{ \dd s^\nu}{\dd t}\right)\right) - \hat{P} \exp\left(ig\int_0^1\dd t A_\mu(s) \frac{ \dd s^\mu}{\dd t}\right).
\end{align}
In order to simplify what follows, let's define $F(t) \equiv igA_\mu(s)\frac{ \dd s^\mu}{\dd t}$ and $G(t) \equiv ig\frac{ \dd}{\dd t} (A_\mu(s) \delta s^\mu) + igF_{\mu\nu} \delta s^\mu \frac{ \dd s^\nu}{\dd t}$. Then we have that $F(t)$ is $\mathcal{O}(\delta^0)$ and $G(t)$ is $\mathcal{O}(\delta)$. 
\begin{align} \label{deltav}
\delta V = \hat{P} \exp\left( \int_0^1 \dd t\left( F(t) + G(t)\right)\right) - \hat{P} \exp\left(\int_0^1\dd t \,F(t)\right).
\end{align}
Expanding the first exponential, we then have
\begin{align}
\hat{P}\left( 1+ \int \dd t(F(t) + G(t)) + \frac12\int \dd t^\prime \dd t(F(t) + G(t))(F(t^\prime) + G(t^\prime)) + ...\right) 
\end{align}
Then let's cancel any terms of order $\delta^2$ or higher. This means that we can only keep terms with $F(t)'s$ only and terms with only a single $G(t)$. Due to the path ordering operator, we can treat all the matrices inside the parentheses as commuting. Thus, it simplifies to:
\begin{align}
\hat{P}\left( \exp\left( \int \dd t F(t)\right)  + \int \dd t G(t)\left( 1+ \int \dd t^\prime F(t^\prime) + \frac12\int \dd t^\prime \dd t^{\prime\prime} F(t^{\prime})F(t^{\prime\prime}) + ...\right)\right).
\end{align}
Putting this all back as the first term in Eq. \ref{deltav}, we see that we're left with:
\begin{align}
\delta V = \hat{P} \int \dd t\, G(t) \exp\left( \int \dd t^\prime\, F(t^\prime)\right).
\end{align}
Now we can put back in the definition of $F(t)$ and $G(t)$:
\begin{align}
\delta V = \hat{P} \int \dd t\, ig\left( \frac{ \dd}{\dd t} (A_\mu(s)\delta s^\mu) + F_{\mu \nu} \delta s^\mu \frac{ \dd s^\nu}{\dd t}\right)\exp\left(ig \int \dd t^\prime A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime}\right).
\end{align}
Notice that the exponential at the end here is just the Wilson line itself, except that now it needs to be path ordered along with the first integral over $t$. Let's now split the first integral into its two terms. The first integral is simply 
\begin{align}
\int \dd t \frac{ \dd}{\dd t}(A_\mu(s)\delta s^\mu) = A_\mu(s(1))\delta s^\mu(1) - A_\mu(s(0))\delta s^\mu(0).
\end{align}
Which can be rewritten in another integral form:
\begin{align}
\int \dd t\, A_\mu(s)\delta s^\mu ( \delta(t-1) - \delta(t)).
\end{align}
We are being sloppy with the bounds on the integrations, but we could always go back to the beginning and define the integral to be over the range $t \in (-\infty,\infty)$ with appropriate step functions to fix the bounds. Now we have
\begin{align}
\delta V = \hat{P}ig \int \dd t\left(A_\mu(s) (\delta(t-1) - \delta(t)) + F_{\mu \nu} \frac{ \dd s^\nu}{\dd t}\right) \exp\left(ig\int \dd t^\prime \, A_\sigma(s)\frac{ \dd s^\sigma}{\dd t^\prime}\right) \delta s^\mu.
\end{align}
Comparing with Eq. \ref{funcderiv} (Sorry, I'm not sure how to get hyperref to work here), we can now read off the functional derivative.
\begin{align}
\frac{ \delta V}{\delta s^\mu (t)}  = \hat{P}\,ig\left[ A_\mu(s)(\delta(t-1) - \delta(t)) + F_{\mu \nu} \frac{ \dd s^\nu}{\dd t}\right] \exp\left(ig\int \dd t^\prime \, A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime} \right).
\end{align}
We can now insert this result into Eq. \ref{Dform}(sorry again!) to get the partial derivative of the Wilson line:
\begin{align}
\frac{ \D V(z,y)}{\D x^\lambda} = \hat{P}\,ig\int \dd t \left[ A_\mu(s)(\delta(t-1) - \delta(t)) + F_{\mu \nu} \frac{ \dd s^\nu}{\dd t}\right] \exp\left(ig\int \dd t^\prime \, A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime} \right) \frac{ \D s^\mu}{\D x^\lambda}
\end{align}
Here we can make some simplifications. First, we can separate the terms and perform the $t$ integral in two of the terms using the delta functions.
\begin{align}
\frac{ \D V(z,y)}{\D x^\lambda} =& \hat{P}\,ig A_\mu(z)\frac{ \D z^\mu}{\D x^\lambda}\exp\left( \int \dd t^\prime A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime}\right) - \hat{P}\,ig A_\mu(y) \frac{ \D y^\mu}{\D x^\lambda} \exp\left( \int \dd t^\prime A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime}\right)\\ &+\hat{P}\,ig \int \dd t\, F_{\mu \nu} \frac{ \dd s^\nu}{\dd t} \frac{ \D s^\mu}{\D x^\lambda} \exp\left( ig\int \dd t^\prime A_\sigma(s) \frac{ \dd s^\sigma}{\dd t^\prime}\right).
\end{align}
Since $z$ is already the end point of the path, in the first term we can pass the path ordering on to the exponential, thus turning it into the Wilson line itself. In the second term, since $y$ is the starting point, we move $A_\mu(y)$ to the right of the exponential, and then put the path ordering operator on the exponential alone, thus also creating the Wilson line. In the third term, let's "un-parametrize" the path so that $s(t) \rightarrow s$ and $s(t^\prime) \rightarrow s^\prime$. 
\begin{align}
\frac{ \D V(z,y)}{\D x^\lambda} =ig A_\mu(z) V(z,y) \frac{ \D z^\mu}{\D x^\lambda} -ig V(z,y)A_\mu(y)\frac{ \D y^\mu}{\D x^\lambda} + \hat{P}\,ig\int_y^z F_{\mu \nu}(s) \dd s^\nu \exp\left(ig \int_y^z A_\sigma(s^\prime) \dd s^{\sigma^\prime}\right)\frac{ \D s^\mu}{\D x^\lambda}.
\end{align}
In the third term, we can account for the path ordering by splitting the exponential into two pieces, the part of the path before $s$, and the part of the path after $s$. So,
\begin{align}
&\hat{P} \,ig\int_y^z F_{\mu \nu}(s) \dd s^\nu \exp\left(ig \int_y^z A_\sigma(s^\prime) \dd s^{\sigma^\prime}\right) \frac{ \D s^\mu}{\D x^\lambda}  \\ &= ig\int_y^z \dd s^\nu \left[ \hat{P}\exp\left(ig\int_s^z A_\sigma(s^\prime) \dd s^{\sigma^\prime}\right)\right] F_{\mu \nu}(s) \frac{ \D s^\mu}{\D x^\lambda} \left[ \hat{P}\exp\left(ig\int_y^s A_\tau(s^\prime)\dd s^{\tau^\prime}\right) \right]\\
&= ig\int_y^z \dd s^\nu\, V(z,s)F_{\mu \nu}(s) \frac{ \D s^\mu}{\D x^\lambda} V(s,y).
\end{align}
We can also write the first two terms in a simplified way:
\begin{align}
igA_\mu(z) V(z,y) \frac{ \D z^\mu}{\D x^\lambda} - igV(z,y)A_\mu(y)\frac{ \D y^\mu}{\D x^\lambda}  = igV(z,s)A_\mu(s) \frac{ \D s^\mu}{\D x^\lambda} V(s,y)|^{s=z}_{s=y}.
\end{align}
Where we've exploited the fact that $V(z,z) = V(y,y) = 1$. Thus our final result is as follows.
\begin{align}
\frac{ \D V(z,y)}{\D x^\lambda} =  igV(z,s)A_\mu(s) \frac{ \D s^\mu}{\D x^\lambda} V(s,y)|^{s=z}_{s=y} +  ig\int_y^z \dd s^\nu\, V(z,s)F_{\mu \nu}(s) \frac{ \D s^\mu}{\D x^\lambda} V(s,y).
\end{align}
Where the first term is used if $x$ is either the start point or the end point for the Wilson line, and the last term is used if $x$ is an interior point of the Wilson line. 
