diff --git a/src/eko/scale_variations/expanded.py b/src/eko/scale_variations/expanded.py
index aeea77b76..7118b03fa 100644
--- a/src/eko/scale_variations/expanded.py
+++ b/src/eko/scale_variations/expanded.py
@@ -1,5 +1,5 @@
 # -*- coding: utf-8 -*-
-r"""This module contains the scale variation operator for the expanded scheme (``ModSV=expanded``).
+r"""Contains the scale variation operator for the expanded scheme (``ModSV=expanded``).
 
 The expressions can be obtained using Eqs. (3.33) and (3.38) of :cite:`AbdulKhalek:2019ihb`.
 Be aware that corresponding the signs of the ingredients there are a number of differences.
@@ -15,7 +15,7 @@
 
 @nb.njit(cache=True)
 def variation_as1(gamma, L):
-    r"""Computes the |NLO| anomalous dimension variation.
+    r"""Compute the |NLO| anomalous dimension variation.
 
     Parameters
     ----------
@@ -34,7 +34,7 @@ def variation_as1(gamma, L):
 
 @nb.njit(cache=True)
 def variation_as2(gamma, L, beta0, g0e2):
-    r"""Computes the |NNLO| anomalous dimension variation.
+    r"""Compute the |NNLO| anomalous dimension variation.
 
     Parameters
     ----------
@@ -57,7 +57,7 @@ def variation_as2(gamma, L, beta0, g0e2):
 
 @nb.njit(cache=True)
 def variation_as3(gamma, L, beta0, beta1, g0e2, g0e3, g1g0, g0g1):
-    r"""Computes the |N3LO| anomalous dimension variation.
+    r"""Compute the |N3LO| anomalous dimension variation.
 
     Parameters
     ----------
@@ -118,14 +118,14 @@ def non_singlet_variation(gamma, a_s, order, nf, L):
     """
     sv_ker = 1.0
     if order[0] >= 2:
-        sv_ker += a_s * variation_as1(gamma, L)
+        sv_ker -= a_s * variation_as1(gamma, L)
     if order[0] >= 3:
         beta0 = beta.beta_qcd_as2(nf)
-        sv_ker += a_s**2 * variation_as2(gamma, L, beta0, gamma[0] ** 2)
+        sv_ker -= a_s**2 * variation_as2(gamma, L, beta0, gamma[0] ** 2)
     if order[0] >= 4:
         beta1 = beta.beta_qcd((3, 0), nf)
         g0g1 = gamma[0] * gamma[1]
-        sv_ker += a_s**3 * variation_as3(
+        sv_ker -= a_s**3 * variation_as3(
             gamma, L, beta0, beta1, gamma[0] ** 2, gamma[0] ** 3, g0g1, g0g1
         )
     return sv_ker
@@ -156,18 +156,18 @@ def singlet_variation(gamma, a_s, order, nf, L):
     sv_ker = np.eye(2, dtype=np.complex_)
     gamma = np.ascontiguousarray(gamma)
     if order[0] >= 2:
-        sv_ker += a_s * variation_as1(gamma, L)
+        sv_ker -= a_s * variation_as1(gamma, L)
     if order[0] >= 3:
         beta0 = beta.beta_qcd_as2(nf)
         gamma0e2 = gamma[0] @ gamma[0]
-        sv_ker += a_s**2 * variation_as2(gamma, L, beta0, gamma0e2)
+        sv_ker -= a_s**2 * variation_as2(gamma, L, beta0, gamma0e2)
     if order[0] >= 4:
         beta1 = beta.beta_qcd((3, 0), nf)
         gamma0e3 = gamma0e2 @ gamma[0]
         # here the product is not commutative
         g1g0 = gamma[1] @ gamma[0]
         g0g1 = gamma[0] @ gamma[1]
-        sv_ker += a_s**3 * variation_as3(
+        sv_ker -= a_s**3 * variation_as3(
             gamma, L, beta0, beta1, gamma0e2, gamma0e3, g1g0, g0g1
         )
     return sv_ker
diff --git a/tests/eko/test_sv_expanded.py b/tests/eko/test_sv_expanded.py
index eae6e9b96..0f49f0442 100644
--- a/tests/eko/test_sv_expanded.py
+++ b/tests/eko/test_sv_expanded.py
@@ -91,7 +91,12 @@ def scheme_diff(g, k, pto, is_singlet):
             )
             ker_b = ker * expanded.non_singlet_variation(gns, a1, order, nf, L)
             ns_diff = scheme_diff(gns, L, order, False)
-            np.testing.assert_allclose((ker_a - ker_b) / ker, ns_diff, atol=1e-3)
+            np.testing.assert_allclose(
+                (ker_a - ker_b) / ker,
+                ns_diff,
+                atol=1e-3,
+                err_msg=f"L={L},order={order},non-singlet",
+            )
 
             # Singlet kernels
             gs = gamma_singlet(order, n, nf)
@@ -105,5 +110,8 @@ def scheme_diff(g, k, pto, is_singlet):
             ker_b = ker @ expanded.singlet_variation(gs, a1, order, nf, L)
             s_diff = scheme_diff(gs, L, order, True)
             np.testing.assert_allclose(
-                (ker_a - ker_b) @ np.linalg.inv(ker), s_diff, atol=5e-3
+                (ker_a - ker_b) @ np.linalg.inv(ker),
+                s_diff,
+                atol=5e-3,
+                err_msg=f"L={L},order={order},singlet",
             )