In theories of gravity, the Riemannian/metric connection, is allowed to have torsion, of which the Levi-Civita connection is the particular torsion-free case.
In the gauge theoretic description of QFTs, the gauge connection defines the covariant derivative and the curvature associated to this connection is the field strength.
Is there also a notion of torsion associated to the gauge connection?
There isn't.
Torsion represents the failure of "infinitesimal parallelograms" to close.
To see this, consider a point $x$ and and two "infintesimally nearby" points $x_1^\mu=x^\mu+\frac{dx^\mu}{d\tau}d\tau=x^\mu+T^\mu d\tau$ and $x_2^\mu =x^\mu+S^\mu d\sigma$.
Now we parallel transport $T^\mu d\tau$ to the second point $x^\mu+S^\mu d\sigma$:
$$ T^\mu(x_2)d\tau=(T^\mu-\Gamma^\mu_{\nu\sigma}S^\nu d\sigma T^\sigma) d\tau, $$ add this to the second point and get $$x_3^\mu= x^\mu +S^\mu d\sigma+T^\mu d\tau-\Gamma^\mu_{\nu\sigma}S^\nu T^\sigma d\sigma d\tau. $$
What if we do the other way around and transport $S^\mu d\sigma$ to $x_1$?
We get $$ \tilde{x}_3^\mu=x^\mu+T^\mu d\tau+S^\mu d\sigma-\Gamma^\mu_{\nu\sigma}T^\nu S^\sigma d\sigma d\tau. $$
The difference of these two points is $$\Delta x_3^\mu=(\Gamma^\mu_{\nu\sigma}-\Gamma^\mu_{\sigma\nu})d\sigma d\tau=T^\mu_{\ \nu\sigma}d\sigma d\tau.$$
The reason I say that is to illustrate that the core concept of torsion involves parallel transporting one vector along another, and then doing the same operation in reverse, transporting the "other" vector along the former one. The invariant definition $$ T(X,Y)=\nabla_XY-\nabla_YX-[X,Y] $$also implies this.
For a general connection, the direction of parallel transport is given by a vector, but the object to be parallel transported are not vectors but general fiber elements. So the operation of transporting a fiber element along a vector, and then transporting the vector along a fiber element makes no sense, because a fiber element doesn't determine a direction in the manifold in general. Eg. the expression $D_Xs$ cannot be symmetrized or antisymmetrized, since $X$ and $s$ are very different objects.
Curvature doesn't have the same problem, since there two vectors are needed to generate a loop, but only a fiber element gets transported.
No, because the torsion of a connection is defined as the action of this connection on the solder form.
$$T = D \theta$$
The Yang-Mill connection has no solder form associated to it.
The generalization is that the symmetry underlying the gauge fields acts on other fields and the gauge-covariant derivatives of those other fields are analogous to the torsion. An example is the Higgs field in the Standard Model.
First: some conventions that will make it a lot easier to provide more detail in the answer to the question.
For the following, I'll use indexed notation for linear functionals, e.g. $T^α_{UV}$ is a functional of a co-vector $α = α_ρ dx^μ$ and vectors $U^μ ∂_μ$ and $V^ν ∂_ν$: co-vectors go upstairs and vectors downstairs. The usual index notation can be fit into this in an upwardly-compatible way by, e.g. $T^ρ_{μν} = T^{dx^ρ}_{∂_μ ∂_ν}$. The functional is "tensorial" in an index if it is linear in that index, and it is a tensor if it is linear in all of its indices. The Kronecker delta becomes the contraction operator $δ^α_U = δ^ρ_μ α_ρ U^μ = α_ν U^ν = α·U$, and is tensor. The connection comes out of the covariant derivative operator by $$Γ^α_{UV} = \left(∇_U V\right)·α = α_ρ U^μ \left(δ^ρ_ν ∂_μ + Γ^ρ_{μν}\right) V^ν$$ and is only tensorial in $α$ and $U$, but not in $V$. As another example, the "structure coefficients" functional (for lack of an official name) $$f^α_{UV} = [U,V]·α = α_ρ \left(U^μ ∂_μ V^ρ - V^ν ∂_ν U^ρ\right)$$ is only tensorial in $α$. For coordinate indices, it vanishes: $f^ρ_{μν} = 0$. An example of an operator-valued functional that is tensorial is the directional derivative operator: $$∂_U = U^ρ ∂_ρ,$$ which is treated, conventionally, in the mathematical literature as the corresponding co-vector, itself: $∂_U = U$.
Now, start from the Cartan structure equations $$Ω^a_b = dω^a_b + ω^a_c ∧ ω^c_b,\quad Θ^a = dθ^a + ω^a_c ∧ θ^c.$$ The players involved in this are the following.
The connection and curvature are associated with the gauge group $GL(n)$, where $n$ is the dimension of the underlying manifold. They can be combined into the respective Lie-valued forms: $$A = ω^c_a E^a_c,\quad F = Ω^c_a E^a_c$$ with the following Lie brackets $$\left[E^a_c, E^b_d\right] = δ^a_d E^b_c - δ^b_c E^a_d.$$ In component form, the one-form and two-form are both Lie-valued: $$A = A_μ dx^μ\quad (A_μ = Γ^c_{μa} E^c_a),\quad F = ½ F_{μν} dx^μ ∧ dx^ν\quad (F_{μν} = R^c_{aμν} E^c_a),$$ with $$F_{μν} = ∂_μ A_ν - ∂_ν A_μ + \left[A_μ, A_ν\right].$$
The frame and torsion fields can be combined into vector-valued one-forms and two-forms: $$q = θ^a E_a,\quad v = Θ^a E_a,$$ with vector-valued component forms $$q = q_μ dx^μ\quad (q_μ = h^a_μ E_a),\quad v = ½ v_{μν} dx^μ ∧ dx^ν\quad (v_{μν} = T^a_{μν} E_a).$$ The basis of the vector space is given by the $E_a$, and they are acted on to the left by the Lie algebra, with $$E^c_a · E_b = δ^c_b E_a.$$ Correspondingly, we can write $$v_{μν} = ∂_μ q_ν - ∂_ν q_μ + A_μ · q_ν - A_ν · q_μ.$$
This is what we'll generalize to other gauge groups and for fields $q$ that are form-valued.
To generalize it, rewrite everything algebraically by adopting the conventions that the basis vectors $\left(E^a_c\right)$ and $\left(E_c\right)$ freely intersperse with the $dx^μ$, that $\left[E^a_c, E^b_d\right] = E^a_c E^b_d - E^b_d E^a_c$ (i.e. by embedding the Lie algebra of $GL(n)$ in a covering algebra), that $E^a_c E_d = E^a_c · E_d$ and that differential forms multiply by the wedge product. Then, using the anti-symmetry $dx^μ ∧ dx^ν = -dx^ν ∧ dx^μ$ we have $$\begin{align} A^2 &= \left(A_μ dx^μ\right) \left(A_ν dx^ν\right)\\ &= A_μ A_ν dx^μ ∧ dx^ν\\ &= ½ \left(A_μ A_ν dx^μ ∧ dx^ν + A_ν A_μ dx^ν ∧ dx^μ\right)\\ &= ½ \left(A_μ A_ν - A_ν A_μ\right) dx^μ ∧ dx^ν\\ &= ½ \left[A_μ, A_ν\right] dx^μ ∧ dx^ν. \end{align}$$ Therefore $$d A + A^2 = F.$$ The corresponding "Bianchi identity" can then be readily derived: $$\begin{align} d F + A F &= d(d A) + d(A^2) + A (d A) + A (A^2)\\ &= 0 + (d A A - A d A) + A d A + A^2 A\\ &= (d A + A^2) A\\ &= F A. \end{align}$$ Thus $$d F + A F - F A = 0.$$
In addition, we have $$\begin{align} A q &= \left(A_μ dx^μ\right) \left(q_ν dx^ν\right)\\ &= A_μ · q_ν dx^μ ∧ dx^ν\\ &= ½ \left(A_μ · q_ν dx^μ ∧ dx^ν + A_ν · q_μ dx^ν ∧ dx^μ\right)\\ &= ½ \left(A_μ · q_ν - A_ν · q_μ\right) dx^μ ∧ dx^ν, \end{align}$$ making use of the anti-symmetry of the wedge product, again. Thus $$d q + A q = v.$$
The analogue of the Bianchi identity is given by $$\begin{align} d v + A v &= d (d q) + d (A q) + A (d q) + A (A q)\\ &= 0 + (d A q - A d q) + A d q + A^2 q\\ &= (d A + A^2) q\\ &= F q. \end{align}$$
As an example, for a gauge field $A = A^a_μ Y_a dx^μ$, with a Lie basis $\left(Y_a\right)$ and Lie brackets $\left[Y_a, Y_b\right] = f^c_{ab} Y_c$, the field strength $F = ½ F^c_{μν} dx^μ ∧ dx^ν$ is given component-wise by $$F^c_{μν} = ∂_μ A^c_ν - ∂_ν A^c_μ + f^c_{ab} A^a_μ A^b_ν.$$ Suppose there are another set of fields - scalar fields - $q = q^A e_A$ that are acted on the left by $Y_a · e_B = f^C_{aB} e_C$. The two sets of indices $(a)$ and $(A)$ may need not have the same dimension.
Then $v = v^C_μ e_C dx^μ$, where $$v^C_μ = ∂_μ q^C + f^C_{aB} A^a_μ q^B.$$ This $v$ field is the gauge-covariant derivative of the $q$ field.
As a more concrete case-in-point: the gauge field $A$ is the boson field for the Standard Model, the gauge group is $U(1) × SU(2) × SU(3)$, and the scalar field $q$ is the Higgs.
