1) Let us write the Wilson-line of a [simple open curve](https://en.wikipedia.org/wiki/Simple_curve#Topology) $\gamma: [s_i,s_f]\to \mathbb{R}^4$ as
$$ U(s_f,s_i)
~=~ \mathcal{P}\exp \left[ i\int_{\gamma} A_{\mu}~ dx^{\mu} \right].\tag{1} $$

2) The [path-ordering](https://en.wikipedia.org/wiki/Path-ordering) $\mathcal{P}$ becomes important if the gauge potential
$$A_{\mu}~=~A^a_{\mu} T_a\tag{2}$$
is non-abelian. Here $T_a$ are the generators of the corresponding Lie algebra.

3) The Wilson-line has [groupoid](https://en.wikipedia.org/wiki/Groupoid) properties, e.g.,
$$U(s_3,s_2)U(s_2,s_1)~=~ U(s_3,s_1), \qquad  U(s,s) ~=~ {\bf 1}.\tag{3}$$

4) If one differentiates wrt. the final point $s_f$, one gets
$$\frac {\mathrm{d}U(s_f,s_i)}{\mathrm{d}s_f} 
~=~  i\dot{\gamma}^{\mu}(s_f)~A_{\mu}(\gamma(s_f)) ~U(s_f,s_i). \tag{4}$$

5) If one differentiates wrt. the initial point $s_i$, one gets
$$ \frac {\mathrm{d}U(s_f,s_i)}{\mathrm{d}s_i} 
~=~  -U(s_f,s_i)~i\dot{\gamma}^{\mu}(s_i)~A_{\mu}(\gamma(s_i)) . \tag{5}$$

6) OP wants to differentiate the Wilson-line $U(s_f,s_i)$ functionally wrt. the gauge potential components $A^a_{\mu}(x)$. One gets
$$ \frac {\delta U(s_f,s_i)}{\delta A^a_{\mu}(x)} 
~=~\int_{s_i}^{s_f}\! \mathrm{d}s~ U(s_f,s)~ 
i\dot{\gamma}^{\mu}(s)\delta^4(x-\gamma(s))T_a~U(s,s_i).
\tag{6}\label{eq:6}$$

7) _Heuristic proof of $\eqref{eq:6}$._ Since we have already used the letter $x\in\mathbb{R}^4$ in $\eqref{eq:6}$ as a fixed space-time point, let us call an arbitrary spacetime point for $y\in\mathbb{R}^4$. 

   - Imagine that $\tilde{A}(y)=A(y)+\delta A(y)$ is an infinitesimal variation of the gauge potential $A(y)$.

   - Imagine that $\delta A(y)$ only differs from zero in an infinitesimally small neighborhood $\Omega$ of the fixed space-time point $x$. 

   - Assume that the curve $\gamma$ intersects the neighborhood $\Omega$ at the parametervalue interval $[s_x-\varepsilon,s_x+\varepsilon]\subseteq [s_i,s_f]$. (If the curve $\gamma$ does not intersect the neighborhood $\Omega$, then the equation $\eqref{eq:6}$ becomes trivially correct: $0=0$.)

   On one hand, such infinitesimal variation of the gauge potential yields
$$\delta U(s_f,s_i)~=~U(s_f,s_x+\varepsilon)~\delta U(s_x+\varepsilon,s_x-\varepsilon)~U(s_x-\varepsilon,s_i), \tag{7}\label{eq:7}$$
and 
$$\begin{align}\delta U(s_x+\varepsilon,s_x-\varepsilon)~\approx~&2\varepsilon i~ \dot{\gamma}^{\mu}(s_x)~\delta A_{\mu}(\gamma(s_x)) \cr
~=~&\int_{\Omega} \!\mathrm{d}^4y~\delta^4(y-\gamma(s_x))~2\varepsilon i\dot{\gamma}^{\mu}(s_x)~\delta A_{\mu}(y)\cr
~\approx~& \int_{\Omega} \!\mathrm{d}^4y~\int_{s_x-\varepsilon}^{s_x+\varepsilon}\! \mathrm{d}s~\delta^4(y-\gamma(s))~i\dot{\gamma}^{\mu}(s)~\delta A_{\mu}(y).\end{align}\tag{8}\label{eq:8}$$ 
On the other hand, the defining property of a functional derivative yields
$$\begin{align}\delta U(s_f,s_i) ~=~&\int_{\mathbb{R}^4} \!\mathrm{d}^4y~ \frac {\delta U(s_f,s_i)}{\delta A^a_{\mu}(y)} ~\delta A^a_{\mu}(y)\cr 
&~=~\int_{\Omega} \!\mathrm{d}^4y~ \frac {\delta U(s_f,s_i)}{\delta A^a_{\mu}(y)} ~\delta A^a_{\mu}(y).\end{align}\tag{9}\label{eq:9}$$
A comparison of eqs. $\eqref{eq:7}$, $\eqref{eq:8}$ and $\eqref{eq:9}$ yields eq. $\eqref{eq:6}$.