Integrating Out New Fermions at One Loop

We present the fermionic universal one--loop effective action obtained by integrating out heavy vector--like fermions at one loop using functional techniques. Even though previous approaches are able to handle integrating out heavy fermions with non--chiral interactions, i.e. vanishing $\gamma^5$ interaction terms, the computations proceed in a tedious manner that obscures a physical interpretation. We show how directly tackling the fermionic functional determinant not only allows for a much simpler and transparent computation, but is also able to account for chiral interaction terms in a simple, algorithmic way. Finally, we apply the obtained results to integrate out at one loop the vector--like fermions appearing in a toy model and in a fermionic model that exhibits strong cosmological phase transitions.


Introduction
Despite a huge experimental effort, signals of physics beyond the Standard Model (SM) remain elusive up to date in direct searches at the Large Hadron Collider (LHC). The absence of such signals has determined a shift of attention towards possible indirect effects of heavy new particles, which can be systematically studied in the context of effective field theories (EFTs) such as the Standard Model Effective Theory (SMEFT) [1,2].
The main strength of EFTs is model independence: experimental measurements can be used to place bounds on the higher-dimensional deformations of a given low-energy theory (e.g. the SM) without specifying any underlying ultraviolet (UV) theory responsible for inducing those deformations. Conversely, the heavy degrees of freedom of a given UV theory can be integrated out and matched to a low-energy EFT, allowing for an efficient study of the resulting effects. In this sense, EFTs serve as a bridge between UV theories and experiments, as long as the new degrees of freedom are sufficiently heavier than the energy scale of the experiment.
Given the current (and future) experimental precision, an accurate translation of the SMEFT bounds to bounds on specific new physics scenarios requires the matching to be performed at the level of (at least) one loop. Traditionally, this task is done at the amplitude level with the help of Feynman (loop) diagrams. A perhaps more elegant and simpler alternative for performing the one-loop matching relies on working directly with the path integral. The idea behind this approach is to identify the contributions of the heavy fields to the one-loop functional determinant, and then expand the determinant in inverse powers of the heavy masses to obtain effective operators containing the light fields. Among the desirable features of such functional methods, there are at least two worth mentioning. Firstly, unlike in the case of Feynman diagrams, gauge covariance is preserved in the intermediate steps by performing a covariant derivative expansion (CDE), which automatically ensures a gauge-invariant final result. Secondly, such methods give universal results that can be applied to a broad class of new physics scenarios with almost no assumption regarding the UV dynamics.
Functional techinques for one-loop matching and their applications were first developed almost 40 years ago in Refs. [3][4][5]. The subject has been recently brought back into focus by Ref. [6], which provided a universal master formula for one-loop matching assuming degenerate masses for the new particles. The generalization to non-degenerate spectra was completed in Ref. [7], and the resulting master formula was named the "Universal One-Loop Effective Action" (UOLEA). However, as pointed as out in Ref. [8] (see also Ref. [9]), both master formulas presented in Refs. [6,7] could only account for the so-called "heavy-only" terms, i.e. terms stemming from loops containing only heavy fields. Shortly afterwards, master formulas including also the "(mixed) heavy-light" terms (originating from loops containing both light and heavy fields) were computed using various methods in Refs. [10][11][12]. With the help of the covariant diagram technique developed in Ref. [13], Ref. [14] established that heavy-only and heavy-light terms share the same structure, and extended the UOLEA of Ref. [7] to include heavy-light terms as well. In addition to the references mentioned up to now, applications of functional methods for one-loop matching have been studied e.g. in Refs. [15][16][17][18][19][20][21] (see also Refs. [22,23]).
Although some steps have been already taken in this direction in Ref. [24], the universal master formulas available to date in the literature do not systematically capture the effects of integrating out heavy fermions at one loop. The main reason behind this is the presence of fermionic interaction terms containing the γ 5 Dirac matrix. Furthermore, even if the terms containing γ 5 are set to zero, applying the existing master formulas to the fermionic case proves to be a tedious task.
It is therefore the aim of this paper to provide a universal master formula that is suitable for integrating out heavy fermions at one loop. To this end, we calculate for the first time the heavy-only contributions arising from integrating out heavy vector-like fermions (VLFs) at one-loop, consistently taking into account the effects of γ 5 . We choose to perform this operation in the unbroken (symmetric) phase, where none of the light scalars have a vacuum expectation value (VEV). In the familiar case where the low-energy theory is the SM, this amounts to calculating the SMEFT dimension-6 operators and Wilson coefficients induced by loops containing only heavy VLFs.
Our work is outlined as follows. In Sec. 2, we specify the general Lagrangian containing all the VLF interaction terms allowed in the unbroken phase, and then set up the expansion of the functional determinant in inverse powers of the VLF masses. We also provide a brief review of alternative methods of tackling the fermionic functional determinant. Sec. 3 is dedicated to the calculation of the fermionic UOLEA terms and their associated universal coefficients. The results from this section are connected to previous universal results in App. C, and then summarized in App. D. In Sec. 4, we apply the results obtained in Sec. 3 to a fermionic toy model and a more realistic VLF model. Finally, we conclude in Sec. 5.

Setup
We consider a fermionic model containing several heavy Dirac VLF multiplets Ψ i , charged under a generic (semi-simple) gauge group, e.g. the SM gauge group SU(3) C ⊗ SU(2) L ⊗ U(1) Y . Working in the unbroken phase of the theory, the most general renormalizable and gauge-invariant Lagrangian of the VLF sector reads where summing over the multiplet (flavour) indices i, j is implicit. Throughout this paper we use the latin indices {i, j, k, l, m, n} to denote flavour indices and assume that summing over them is implicit when considering Lagrangians. 1 By virtue of working in the unbroken phase, the vector-like mass and heavy fermion covariant derivative matrices are diagonal in multiplet space, M ij = m i δ ij and D µ ij = D µ i δ ij (no sum), whereas the S and P matrices are in general non-diagonal, which implies that they do not always commute with M. Hermicity of the Lagrangian in Eq. (1) implies that both S and P are Hermitian.
The interactions of the heavy fermion fields Ψ i with the light gauge fields are encoded in the D µ matrix, which takes a generic form where the sum runs over all the light gauge fields G, and the repeated index a implies summing over all generators of each gauge group. The index i of the group generator t a R,i serves as a reminder that it lies in the same gauge group representation as the heavy fermion Ψ i that it acts upon. On the other hand, S = S(φ) and P = P (φ) specify the Yukawa interactions of Ψ i with light (pseudo)scalar fields denoted generically by φ.
Since we are working in the unbroken phase, where none of the scalars appearing in S and P have a vacuum expectation value (VEV), the vector-like nature of the new fermions Ψ forbids an axial term Ψ i γ µ γ 5 A µ ij Ψ j , while we discard the tensor term Ψ i σ µν T µν ij Ψ j , as it encodes non-renormalizable interaction terms of at least dimension-5. The pseudoscalar term iγ 5 P ij is however allowed, as the left and right chiralities of a given VLF can couple differently to the scalar fields present in S and P . If we were to work in the broken phase, the general Lagragian would contain both the axial term mentiond above and a non-diagonal vector term, Ψ i γ µ V µ ij Ψ j . To obtain the effective action induced by integrating out the heavy fermions Ψ i at one loop, we must calculate the "tracelog" [3][4][5][6] of the operator appearing in Eq. (1): Here, the fermionic factor c f = −1 appears from performing a gaussian path integral over anti-commuting fermionic fields. The symbol Tr denotes a full trace over coordinate space, flavour space, and internal degrees of freedom (spin and gauge). Computing the effective action up to dimension-n amounts to expanding S H up to terms of O(M 4−n ). To set up this expansion, the usual procedure is to write explicitly the trace over coordinate space using momentum eigenstates in d dimensions: where for the last equality we used the fact that "sandwiching" an operator function f (∂ µ ) between the exponential amounts to shifting f (∂ µ ) → f (∂ µ + ip µ ). At this stage, one can "sandwich" the log in Eq. (4) by exp ±iD µ ∂ ∂pµ , as originally done in Refs. [3,5], and obtain an expansion that is manifestly gauge-invariant at all intermediate stages [6], i.e. covariant derivatives appear only in commutators. We do not pursue this avenue in this paper and instead work directly with the expression from Eq. (4).
The trace operator tr in Eq. (4) denotes tracing over flavour, spin, and gauge degrees of freedom. After flipping the sign of the loop momentum, we get: Passing to the second line of Eq. (5) is possible only because the of the equality tr log(A + B) = tr log(A) + tr log(1 + A −1 B). Moreover, we have discarded the constant term ∝ log( / p − M) as it does not depend on any fields, and defined a new momentum integral measure as: Denoting the fermionic propagator as / p − M −1 ≡ / ∆, we can finally perform the Taylor expansion of the log from Eq. (5) and write down the effective one-loop Lagrangian as A desirable feature of the expansion in Eq. (6) is that each order n contains only operators of dimension-n, which renders the power counting transparent. This owes to the fact that all three terms D µ , S, and P are of dimension-1, as they depend linearly on bosonic fields. Therefore, truncating the series at n = 6 ensures the inclusion of all the effective operators arising at dimension-6 or lower.
In the existing literature, several methods have been advanced for evaluating the fermionic functional trace. Refs. [6,7,10,12] use the invariance of the trace under γ µ → −γ µ and Tr log A + Tr log B = Tr log AB to bring the fermionic trace in a form that resembles the bosonic trace: with where we use the usual notation for the (anti)commutator and σ µν = i 2 [γ µ , γ ν ]. Once brought into the form shown in Eq. (7), the fermionic trace can be computed using the known results for the bosonic trace [6,7], dubbed the "Universal One-Loop Effective Action" (UOLEA). However, these results apply solely to the case where U ferm does not contain any open covariant derivatives, 2 which is true only if [γ µ , W ] = 0. Since in the general renormalizable case displayed in Eq. (1) W contains γ 5 , it is clear that the trick shown in Eq. (7) is helpful only if the γ 5 piece is vanishing, i.e. P = 0. And even if P = 0 and therefore [γ µ , W ] = 0, the rather lengthy expression of U ferm from Eq. (8) makes this method impractical for computing the fermionic functional trace.
Another method to compute the fermionic functional trace was put forward in Ref. [13] and relies on a diagrammatic computation of Eq. (6). Exploiting the cyclic property 3 and the gauge invariance of the trace in Eq. (6), one can write down all the allowed operators appearing at a given order and then calculate the corresponding universal coefficients by considering loop diagrams involving combinations of insertions of the three terms from Eq. (6) (iγ µ D µ , S, and iγ 5 P ). This method has been dubbed the "covariant diagram" approach, as the covariant derivative D µ is treated as a single object, as opposed to being split into ∂ µ and the gauge boson piece, as is done in conventional Feynman diagrams. Our computation of the fermionic functional trace will be similar in spirit to the covariant diagram approach, as we also rely on the cyclic property of the trace. However, as opposed to Ref. [13], we do not use diagrams, but read the terms directly from Eq. (6). Moreover, Ref. [13] splits the fermionic propagator into two parts: In our approach, instead of performing this separation, we just compute the fermion traces and the resulting loop integrals using Package-X [25,26].
Finally, we also mention Ref. [24], which discusses the fermionic extension of the UOLEA (both heavy-only and heavy-light contributions, plus mixed scalar-fermion contributions) and calculates all the relevant momentum integrals, but without specifying the γ-matrix structure of the interactions and therefore without evaluating the spin traces.

Computation of the Fermionic One-Loop Effective Action
We now compute order-by-order the terms from the effective lagrangian in Eq. (6), up to n = 6, which corresponds to dimension-6 operators. In the particular case of the SMEFT [27], where the only light scalar field is the SM Higgs doublet, gauge invariance dictates that terms corresponding to odd powers of n vanish. Nevertheless, we focus on the general case and consider odd powers of n as well. Throughout the computation, we denote the trace over spin degrees of freedom (i.e. Dirac traces) as tr s , and the trace over gauge indices as tr g . Moreover, we pull out from each term in the effective action a factor of tr s 1 ≡ n D = 4, which represents the number of spin degrees of freedom for a Dirac fermion (we denote the identity matrix in spinor space as 1). This choice facilitates the comparison with the UOLEA results (see App. C). In addition, we explicitly write down the symmetry factors for each gauge-invariant trace, as opposed to absorbing them in the universal coefficients. For example, the tr g (S ij S jk S ki ) term comes with a factor 1 3 in front, whereas tr g (S ij P jk P ki ) has no symmetry factor. For the universal coefficients, we use the shorthand notation g(m i , m j , m k , m l , · · · ) ≡ g ijkl··· .
Also, we often encounter coefficients with flipped signs for some of the masses. To this end, we define the following notation: i.e. an index between brackets translates to a flipped sign for the corresponding mass. As mentioned before, we consider the sum over flavour indices to be implicit. Finally, before starting our computation, we briefly discuss the type of terms that are allowed in the one-loop effective action. Even if not apparent from Eq. (6), gauge invariance ensures that in the final result all the covariant derivative matrices organize into commutators, 4 i.e. the only dependence on covariant derivatives is through pieces such as [D µ , X] (with X = S, P ) or [D µ , D ν ], which we denote as Moreover, all terms containing an odd number of D µ 's vanish because of Lorentz invariance. Due to the properties of Dirac traces involving γ 5 , the only terms containing odd powers of P are O(P D 4 ) and O(SP D 4 ), which appear at dimension-5 and dimension-6, respectively.

Dimension-1 Terms
At dimension-1, there is only one possible term, O(S) , all the other ones being forbidden by gauge and Lorentz invariance: where g i 1 can be expressed in terms of the master integrals from App. A as The factor n D from Eq. (13) comes from the trivial spin trace tr s 1, and the terms involving γ µ and/or γ 5 matrices vanish under the spin trace. Note that, by gauge invariance, the sole term from Eq. (13) vanishes unless the light particle spectrum contains a real singlet. In this case, this term would represent a tadpole term for the singlet.

Dimension-2 Terms
At dimension-2, only the O(S 2 ) and O(P 2 ) terms appear, as the remaining O(D 2 ) and O(D µ S, D µ P ) terms are forbidden by gauge and Lorentz invariance, respectively. Therefore, we can safely drop the i / D piece from the n = 2 contribution: where we have kept the symmetry factor 1 2 of the gauge traces involved. Both terms quantify the one-loop contribution of the heavy fermions to the masses of the light scalars present in the theory. Using dimensional regularization [28] in d = 4 − 2ǫ dimensions, the universal coefficients are given by: It is interesting to remark that the coefficients of the O(S 2 ) and O(P 2 ) terms are related: g is equal to g ij 2 with the sign of m j flipped. This equivalence follows from the identity (see App. A for the definition ofĝ µν ): where the term proportional toĝ µν is responsible for generating the finite correction δg ij 2 . 5 As detailed in App. A, the term proportional toĝ µν from Eq. (17) stems from using the so-called BMHV scheme to handle γ 5 in d-dimensional space, and it needs to be kept as we are dealing with divergent loop integrals.
Using the equality from Eq. (17) (or variations thereof) and trace symmetry arguments, the coefficients of terms involving even powers of P can be straightforwardly expressed with the help of coefficients of terms involving only S, as illustrated later on. As an added bonus, for terms of dimension 5 and 6 we can setĝ µν → 0, as the corresponding loop integrals are finite, which will greatly simplify our computation.

Dimension-3 Terms
Going forward to n = 3, Lorentz and gauge symmetries restrict the possible terms down to O(S 3 , SP 2 ), which from the physical point of view renormalize the trilinear light scalar couplings present in the unbroken phase. Skipping intermediate steps, the dimension-3 Lagrangian is given by: with the following coefficients: 5 Note that using a scheme in which γ 5 (naively) anticommutes with γ µ in d-dimensional space would imply that δg ij 2 = 0. and g Eq. (11). As pointed out at the end of Sec. 3.2, the coefficient of O(SP 2 ) is related to the coefficient of O(S 3 ) by flipping the sign of m k and adding a finite contribution δg ij 3 , which follows from of Eq. (17). The symmetry factor 1 3 of the gauge trace tr g (S 3 ) comes from its Z 3 symmetry, while the lack of cyclical symmetry of tr g (SP 2 ) explains why it has no symmetry factor.

Dimension-4 Terms
The discussion becomes more involved when passing to dimension-4 or higher terms, as terms with covariant derivatives are now allowed by gauge invariance, unlike the case of n = 1, 2, 3. We organize the three possible terms as where X is used to generically denote S and P . From the physical point of view, the three terms renormalize the quartic scalar couplings, the kinetic terms of the scalars, and the gauge kinetic terms, respectively. As this is the first time we encounter traces containing covariant derivatives, we comment on whether the cyclic property can be used for such traces. When dealing with the trace over internal degrees of freedom only, denoted as tr, the cyclic property obviously does not hold. However, this property does hold when using the full trace, Tr, which includes a trace over the coordinate space in which the derivative operator ∂ µ acts. Since S H = d d x L H , one can convert the trace over internal degrees of freedom to a full trace [4,12] using the identity: where we have used the d-dimensional space-time volume V d to compensate the infinite Dirac distribution, V d = δ d (0) = x|x . Using this trick, one can switch from tr to Tr, apply the cyclic property for covariant derivative terms to cast them into the desired form, and then revert to tr. The net result is that one can safely apply the cyclic property for covariant derivative terms too. This is the reason why it is possible to set up a diagrammatic approach, as done in Ref. [13].
In the present case, however, there are more independent terms 6 involving S and P , namely S 4 , S 2 P 2 , (SP ) 2 , and P 4 . These terms are given by: The universal coefficient g ijkl 4 and the finite corrections δg 4a , δg 4b read while the remaining universal coefficients follow from their definition from Eq. (11). In Eq. (22) the 1 q symmetry factors have been kept in accordance with the Z q symmetry of each trace, while the minus sign in front of the third term is a result of the identity which is a simple variation of Eq. (17). Note that the g term does not receive any finite corrections stemming from the BMHV treatment of γ 5 , as the term proportional tô g µν from Eq. (24) scales as p µ and not p 2 . Consequently,ĝ µν ends up multiplying a finite integral, and the result vanishes when taking the limit ǫ → 0.
O(X 2 D 2 ) terms. For calculating the O(X 2 D 2 ) terms, we use an approach similar to the one presented in Ref. [13]. We focus on the O(S 2 D 2 ) term, from which the remaining O(P 2 D 2 ) term follows immediately, as mentioned at the end of Sec. 3.2.
We first note that there are two independent terms that contain two covariant derivatives and two powers of S, S 2 D 2 and (SD) 2 , and their coefficients follow from the n = 4 term in the Taylor expansion of Eq. (6): Note that each term has the appropriate symmetry factor. Afterwards, we write down all the possible gauge-invariant terms (and the corresponding symmetry factors) arising at O(S 2 D 2 ) -in this case there is only term -and expand the commutators: where we have used the symmetry of the form factor, g ij 5 = g ji 5 , which follows from the symmetry of the associated trace. With Eqs. (25) and (26) at our disposal, we can solve for g ij 5 by equating the factors multiplying the (SD) 2 term: with d = 4 − 2ǫ the number of space-time dimensions. We mention that, equivalently, g ij 5 could have been computed by matching the factors in front of the S 2 D 2 term. While for calculating g ij 5 it is sufficient to equate the prefactors of only one of the two independent terms 7 in Eqs. (25) and (26), considering the other term does have some value, as it provides a consistency check of the results. We have performed this check and re-confirmed Eq. (27).
As for the O(P 2 D 2 ) term, the computation proceeds in an equivalent manner, which allows us to write the effective Lagrangian at O(X 2 D 2 ): where the finite correction is given by O(D 4 ) terms. Since the D µ 's are diagonal in multiplet space, the term involving four covariant derivatives depends only on one mass. The only gauge invariant quantity involving with where the O(ǫ) term in the factor multiplying the divergent integral I[p 4 ] 4 i was retained to correctly account for the finite part. We note that our result for the O(D 4 ) term agrees with the findings of Refs. [6,13].
For clarity, let us break down the intermediate steps in obtaining the result in Eq. (30). Starting from the Taylor expansion in Eq. (6), the O(D 4 ) piece is given by: Passing to the third line is done through a cyclic permutation of the D 2 i D 2 i term, and then Eq. (30) is recovered by noting that

Dimension-5 Terms
For the dimension-5 case, we organize the possible terms as follows: We choose to treat the O(SD 4 ) and O(P D 4 ) terms separately, as they are different in terms of their CP properties.
As pointed out previously, the dimension-5 and dimension-6 universal coefficients are finite, allowing us to compute the corresponding loop integrals in 4 instead of d dimensions. Therefore, γ 5 retains its usual anticommuting properties and we no longer need to keep the pieces proportional toĝ µν from Eqs. (17,24) (or variations thereof) when computing the universal coefficients for terms involving even powers of P . We stress once again that, with the help of Eqs. (17,24) and trace symmetry arguments, these coefficients follow effortlessly from the coefficients of operators containing only insertions of S.
with the universal coefficient given by: terms, we follow the same procedure as for the case of O(X 2 D 2 ) terms. Focusing on the O(S 3 D 2 ) contribution, the only independent gaugeinvariant combination is It is clear from this relation that the easiest way to find g ijk 8 is to compute the loop integral multiplying tr g (S ij S jk D 2 k S ki ) from the covariant derivative expansion in Eq. (6): from which g ijk 8 is found to be 8 Including the terms containing powers of P , the O(X 3 D 2 ) Lagrangian reads: O(SD 4 ) terms. The simplest way to compute this is to directly compute it from the CDE in Eq. (6). The result is: from where the universal coefficient g i 9 can be easily read off as: We chose to write g i 9 explicitly, as expressing it through the master integrals defined in Eq. (A.1) would have lead to a much more cumbersome relation.
Note that the O(SD 4 ) term from Eq. (40) where ε µνρσ is the 4-dimensional Levi-Civita tensor, and we have defined The universal coefficient appearing in Eq. (42) is: which we also write down in terms of master integrals thanks to the compactness of the expression. As in the O(SD 4 ) case, the O(P D 4 ) term vanishes unless there are real (pseudo)scalars present in the light particle spectrum. Contrary to O(SD 4 ), the term in Eq. (42) conserves CP in the case of light pseudoscalars.

Dimension-6 Terms
We now turn our attention towards the final set of terms considered in this work, the dimension-6 terms. Gauge invariance allows for a multitude of n = 6 terms, which we organize as: The O(SP D 4 ) piece is written separately from O(X 2 D 4 ) because it is the only one at this dimension that depends on the dual field strength tensor F µν .
O(X 6 ) terms. Although lengthy, the X 6 piece is straightforward to compute by generalizing from O(X 5 ) and reads: tr g (S ij S jk S kl P lm S mn P ni ) + 1 2 g ijk(l)(m)(n) 11 tr g (S ij S jk P kl S lm S mn P ni ) +g ijk(l)m(n) 11 tr g (S ij S jk P kl P lm P mn P ni ) − g ij(k)(l)m(n) 11 tr g (S ij P jk S kl P lm P mn P ni ) tr g (P ij P jk P kl P lm P mn P ni ) , while the universal coefficient is equal to: To maintain the expression compact, we have defined a ± as: The O(X 6 ) terms suggestively illustrate how the use of trace symmetry and of Eqs. (17) and (24) streamlines our computation. Instead of having to calculate eight operators and their corresponding universal coefficients, it is enough to consider only one operator, S 6 , together its universal coefficient, with the other seven following effortlessly.
O(X 4 D 2 ) terms. As in the case of the other terms with two covariant derivatives, we first compute the universal coefficients for P = 0, and then generalize to P = 0. At O(S 4 D 2 ), we have two gauge invariant terms, which upon expanding the commutators become: In the last equality, we have used the relation g ijkl 13 = g klij 13 (inherited from the symmetry of the associated trace) and performed the symmetrization: To calculate g ijkl 12,13 , we match Eq. (48) on the corresponding terms obtained form the CDE in Eq. (6), which we choose to be Equating the expressions in Eqs. (48) and (49), we find: and As a check, we have also performed the matching on the remaining term, S 2 D µ S 2 D µ , and found it to be consistent with our expressions for g ijkl 12 and g ijkl 13 . Having all the ingredients, we finally write down the full O(X 4 D 2 ) effective Lagrangian, including terms with P as well: Again, trace symmetry arguments together with Eqs. (17, 24) made our task much simpler: instead of fourteen different terms, we only had to compute two, namely S 2 [D µ , S] 2 and (S [D µ , S]) 2 .
O(X 2 D 4 ) terms. Focusing again on the terms containing only S, there are four independent gauge invariant traces that arise at O(S 2 D 4 ). We choose them as follows: where we have expanded the commutators and kept only four open covariant derivative terms, which are necessary for computing the four universal coefficients present at this order. Before matching Eq. (53) to the corresponding terms from the CDE in Eq. (6), we write down some useful relations for calculating the universal coefficients. At O(S 2 D 4 ) (and O(P 2 D 4 )), the loop integrals have four Lorentz indices and thus have the general form 9 a 1 g µν g ρσ + a 2 g µσ g νρ + a 3 g µρ g νσ .
In practice, however, we do not need the full loop integral, but just its scalar components, defined above as a 1,2,3 . To isolate these components, we define P αβδλ = 5g αβ g δλ − g αλ g βδ − g αδ g βλ 72 , such that P µνρσ , P µσνρ , and P µρνσ single out a 1 , a 2 , and a 3 , respectively. The number of space-time dimensions has been set to 4, as all the loop integrals at this order are finite. We now come back to computing the universal coefficients. Comparing Eq. (53) with the relevant terms from Eq. (6), we find: Note that the 1 2 factors in the expressions of g ij 14,16 in Eqs. (54) and (56) are symmetry factors multiplying the SD 2 SD 2 and SD µ D ν SD µ D ν terms in the CDE.
Concerning the terms involving P instead of S, we do not write them explicitly, as they can be obtained in a straightforward manner by taking the gauge-invariant O(S 2 D 4 ) traces from Eq. (53) (first and second lines) and substituting S → P and g ij N → g O(SP D 4 ) terms. At dimension-6, the O(SP D 4 ) terms are the only ones that depend on the dual tensor F µν , therefore we treat them separately. These terms can be computed directly from the CDE in Eq. (6): with the universal coefficients having the rather simple expressions: and F µν defined in Eq. (43). The remaining two terms, SD µ P D ν D ρ D σ and P D µ SD ν D ρ D σ , do not appear because their corresponding loop integrals vanish. This is expected from the physical point of view, as such terms could only originate from an operator such as O(D 6 ) terms. The easiest way of tackling the O(D 6 ) piece is to directly compute the associated loop integral from the effective Lagrangian from Eq. (6), and then group the result in independent terms using the cyclic property of traces: At this order, there are two possible gauge-invariant terms, which we write as: After expanding the commutators from the gauge-invariant traces and comparing Eqs. (60) and (61), we are able to write down the expressions for the two universal coefficients arising at O(D 6 ): We refrain from writing the form factors in terms of master integrals, as the expressions involved would be much lengthier. Moreover, our results for O(D 6 ) are in agreement with the ones previously obtained in Refs. [6,13].

Examples
In this section, we apply our previously obtained results to two concrete scenarios. The first example is a toy model involving a heavy charged fermion coupling to a real pseudoscalar and a (massless) U(1) gauge boson, while the second one represents a fermionic model discussed in Ref. [29] in the context of strongly first order electroweak phase transitions.

A Toy Model
The toy model that we consider contains a massless U(1) gauge boson, A µ , and a light real pseudoscalar, φ, together with a heavy charged VLF ψ, which is to be integrated out at one-loop. Ignoring the pieces that are irrelevant for our purposes, the Lagrangian for this toy theory reads: with λ a real Yukawa coupling, as dictated by the hermicity of the Lagrangian. Assuming unit charge for ψ, the covariant derivative is given by where g is the U(1) gauge coupling. Casting the toy model Lagrangian in the form of Eq. (1), we obtain the following identities: Note that the covariant derivative acting on the real pseudoscalar φ reduces to the usual derivative, as expected. As there is only one heavy fermion to be integrated out, P and [D µ , P ] are scalars in flavour space, which is also true for the field strength tensor, given by: Therefore, the trace in multiplet space is straightforward, and all the universal coefficients depend on only one mass scale, i.e. the heavy fermion mass m. The gauge trace is trivial as well, as both light fields, φ and F µν , have no gauge quantum numbers. We now use the equal mass expressions of the universal coefficients reported in Tables 1-4 to write down the one-loop effective Lagrangian arising from integrating out ψ: It is worthwhile to note the absence of the P [D µ , P ][D ν F νµ ] and F µ ν F ν ρ F ρ µ terms from the results in Eq. (67). The first term vanishes for real (pseudo)scalars such as φ, as can be shown through integration by parts. As for the second term, one can prove that which shows that this effective operator vanishes in the case of abelian gauge fields.

Fermions and Cosmological Phase Transitions
We now focus on applying the techniques from Sec. 3 to a more realistic model in which vector-like (VL) leptons are added to the SM particle spectrum in order to produce a strongly first order EW phase transition [29]. Besides the SM particle content, this model contains three VL lepton 10 multiplets: where we use the notation (SU (3 where P L,R are chiral projectors. The covariant derivatives acting on the fermionic fields read: Although not explicitly written, it is understood that the ∂ µ and B µ pieces in D µ,L are multiplied by the 2 × 2 identity matrix in SU(2) L space. In order to match the VLL Lagrangian in Eq. (70) to the notation used in Eq. (1), we define: with A = E, N. Working in the basis Ψ = L E N T , where T means transposition only in flavour space, the expressions of the S and P matrices read: In the equation above, the subscripts denote the dimension in SU(2) L space of each element of the S and P matrices, while the entries that do not carry any SU(2) L indices have no subscripts. However, from now on, we stop writing the dimensions of each SU(2) L submatrix in order to simplify the notation. Concerning the mass matrix and the covariant derivative matrix, they are given by: To keep the discussion simple, we choose all VLL masses to be degenerate, m L = m E = m N ≡ m, and work in the limit where z E,N = 0, i.e. a vanishing P matrix. 11 Note however that setting P = 0 amounts to removing all CP-odd operators that might arise from integrating out the VLLs. We will return to this subject towards the end of the section. Besides S, the other building blocks appearing in the computation are [D µ , S] and F µν , given by: and Here, the covariant derivatives acting on the Higgs doublet are defined in Eq. (B.4) from App. B, and the SU(2) L and U(1) Y gauge field strengths follow the same notation as in Ref. [27]. As exemplified by Eq. (75), gauge invariance ensures that for any generic matrix X(φ) depending linearly on the light fields φ, we have that [D µ , X(φ)] = X(D µ φ).
Before listing the final result, we provide some details of our calculation. To simplify the computation, it is useful to note that S 2 is a diagonal matrix in flavour space: We now write down the Lagrangian obtained after integrating out the VLLs in Eq. (69) by splitting it into pieces containing dim ≤ 4 operators and another one containing dim = 6 effective operators. The first piece reads: and renormalizes the scalar and gauge kinetic terms, plus the scalar mass and quartic terms. Note that we have used the MS scheme and dropped the divergent parts. Using the SMEFT operator basis defined in App. B, the (CP-even) dimension-6 effective Lagrangian generated by integrating out at one-loop the VLLs in our model is given by: There are several simple consistency checks that one can perform to assess the validity of the results presented in Eq. (80). For example, the leading contribution to the T parameter [30,31] is proportional to the Wilson coefficient of O T : and vanishes in the custodial limit y E = y N , as expected. In the above equation, e is the electromagnetic coupling constant and v the Higgs VEV. We have explicitly checked that our expression for the T parameter from Eq. (81) matches the one obtained in Ref. [32]. Furthermore, the Wilson coefficient of O 6 can alternatively be computed from the Coleman-Weinberg potential [33] corresponding to the fermionic Lagrangian in Eq. (70), and we have explicitly checked that the two methods give the same result. Yet another consistency check can be done by inspecting the physical Higgs boson's loop-induced coupling to photons. The leading O(m −2 ) VLL contribution to the hγγ coupling can be derived from Eq. (80) as well as through low-energy Higgs theorems [34,35] (see also Refs. [36][37][38][39] for VL fermion applications of low-energy Higgs theorems). Denoting the Wilson coefficient of a given operator O X as C X , the VLL contribution to the hγγ coupling can be read from: where e is the electromagnetic coupling constant, A µν the electromagnetic field strength, and Q E = Y − 1 2 is the electrical charge of E (′) . For simplicity, we have set y N = 0 in Eq. (82), as the matching for y N = 0 proceeds in an analogous way. Working in the unitary gauge, the Higgs doublet becomes upon electroweak symmetry breaking, with v the vacuum expectation value and h the physical Higgs scalar.
Using the low-energy theorems (LET), the (CP-even) effective hγγ coupling arising from integrating out the heavy VL leptons from Eq. (69) is described by [39] L LET, CP where M E is the E − E ′ mass matrix and, again, we have set y N = 0 for simplicity. The E-sector mass matrix is extracted from Eq. (70) and reads: Plugging M E in Eq. (83) and expanding up to O(m −2 ), we find the same result as in Eq. (82), which constitutes a useful check for our results in Eq. (80).
Using methods similar to Ref. [6], Ref. [18] performed the same computation as in our Eq. (80). However, our results differ from theirs even after taking into account the redundancies in the operator basis used in Ref. [18]. 12 As advertised earlier, we now turn our attention towards CP-odd dimension-6 operators, which can be generated only if P = 0. Therefore, in the following, we consider z E,N = 0, and focus only on the CP-odd effective Lagrangian. In the VLL model at hand, the CP-odd operators arise from the g ij(k)(l) 13 , g i(j)(k)l 13 terms in Eq. (52) and from the O(SP D 4 ) Lagrangian in Eq. (58). Putting together all these contributions, we find: where again we have set y N L,R = 0. Using the low-energy theorem instead, the VLL contribution to the CP-odd effective hγγ coupling 13 reads: Operators involving gluons would also be generated, the CP-even ones being given by: There is also one CP-odd gluonic operator generated: with O GG defined in App. B.

Summary and Conclusions
Initially proposed more than 30 years ago [3,5], functional methods for integrating out heavy degrees of freedom have been recently revived in Ref. [6] and are currently an ongoing scientific effort [7][8][9][10][11][12][13][14]. Although interesting from the theoretical point of view, the true power of functional methods lies in their practical applications: once a UV sector is specified, matching at one-loop to a low-energy EFT can be achieved by using a few universal master formulas and evaluating matrix traces in internal space (e.g. spin, gauge, flavour).
In this paper, we have extended the universal one-loop formulas presented in Refs. [6,7,[10][11][12][13][14] to include the case of heavy vector-like fermions whose left and right chiralities have different Yukawa interactions with light scalars (e.g. the Higgs boson), as encoded by the S and P matrices in Eq. (1). We have considered the limit where the new fermions do not mix with the SM fermions, the resulting universal coefficients being referred to in the literature as "heavy-only" coefficients. The computations have been performed in the unbroken phase, such that there are no vector and axial current interaction terms for the new fermions (besides the ones encoded in the covariant derivative matrix).
Interestingly, as exemplified throughout Sec. 3, the coefficients of operators containing an even number of P insertions can be easily computed from the corresponding S-only coefficients by flipping the signs of one or more masses and appropriately adjusting the symmetry factors. This led to a drastically simpler computation: out of the 44 universal coefficients arising at dimension-5 and dimension-6, we needed to calculate only 15. The most striking examples are the operators arising at O(X 6 ) and O(X 4 D 2 ). For the former, we had to compute just one coefficient instead of eight, whereas for the latter we had to compute just two coefficients instead of fourteen. As an exception to this rule, in the case of dimension-4 or less coefficients we also had to add a finite correction, which stems from the BMHV treatment of γ 5 matrices in d dimensions.
All in all, we find that the heavy-only fermionic UOLEA computed in the unbroken phase is described by 56 independent operators and their corresponding coefficients, which we summarize in Tables 1-4 in App. D. As a cross-check of our computation, we have computed the S-only coefficients using the universal master formula from Refs. [7,13] and found agreement between the two methods (a "dictionary" between the fermionic universal coefficients and the bosonic UOLEA coefficients is provided in App. C). Finally, we have applied our results in Sec. 4 and integrated out heavy fermions (i) in a toy model with a pseudoscalar Yukawa interaction and (ii) in a more realistic fermionic model which can accommodate a strongly first-order electroweak phase transition.
Concerning future directions, it would certainly be useful to take advantage of the universality of functional methods and derive more general one-loop master formulas that would cover cases involving e.g. open covariant derivatives or mixed statistics (i.e. heavy bosons and fermions integrated out simultaneously). Such developments would considerably simplify phenomenological studies of a broad class of New Physics scenarios, as performing the one-loop matching to e.g. the SMEFT would translate to using a few universal formula(s) and calculating algebraic traces, which also opens up the potential for automation.
work is supported by University of Nebraska-Lincoln, National Science Foundation under grant number PHY-1820891, and the NSF Nebraska EPSCoR under grant number OIA-1557417.
A Master Integrals and Treatment of γ 5 In this appendix, we provide a general definition for the master integrals used to calculate the fermionic one-loop coefficients, and then discuss our treatment of γ 5 in d-dimensional loop integrals. We define the master integrals as: from which one can derive the following relation, which connects our definition to the one in Ref. [13]: with g µ 1 ···µ 2np the completely symmetric tensor, i.e. g µνρσ = g µν g ρσ + g µρ g νσ + g µσ g νρ for n p = 2, and d = 4 − 2ǫ the number of space-time dimensions. Throughout the paper, we use Package-X [25,26] to compute the fermion traces and the resulting master integrals. Analytical expressions of the master integrals can be found in Ref. [13]. We now focus on the problems that arise when dealing with the γ 5 matrix in d = 4 − 2ǫ dimensions. Alongside the Levi-Civita tensor, ε µνρσ , γ 5 is an intrinsically 4-dimensional object and therefore nontrivial to define in d dimensions. In this work, we address this issue by adopting the "Breitenlohner-Maison-'t Hooft-Veltman" (BMHV) scheme [28,41], in which the d-dimensional space is formally separated into a direct sum between a 4dimensional and a (−2ǫ)-dimensional subspace. As a result, each Lorentz vector/tensor now possesses 4-dimensional and (−2ǫ)-dimensional components. Following Ref. [42], we denote the former components by a bar, and the latter by a hat. For example, the ddimensional metric g µν is written as: whereḡ µν andĝ µν reside in 4-dimensional and (−2ǫ)-dimensional spaces, respectively. They act as projectors onto these spaces, i.e.
a µ =ḡ µν a ν ,â µ =ĝ µν a ν ,ā µν =ḡ µρ a ρ ν ,â µν =ĝ µρ a ρ ν , for a generic Lorentz vector/tensor a, and obey the following properties: For a full list of properties, we refer the reader to Ref. [42]. In the following, we just list a few selected relations that are useful for our purposes. We note that relation (A.4) applies to Dirac matrices too, from which we can deduce some of the following properties: where 1 is the identity matrix in spinor space. Finally, γ 5 no longer anticommutes with γ µ , but retains some of its usual 4D properties, such as squaring to identity: As an example, let us calculate the finite correction δg ij 2 appearing in Sec. 3.2. We start by noting that Eq. (A.7) implies from which it follows that where ∆ i,j is defined in Eq. (A.1), and we have omitted the contribution that does not vanish in 4-dimensional space, as it is not relevant for the computation of δg ij 2 . Using the equation above, we calculate δg ij 2 as: where we have replaced p µ p ν → 1 d g µν p 2 under the integral, usedĝ µν g µν =ĝ µνĝ µν = d − 4, and then took the limit ǫ → 0 in the last step.
while the definitions of the field strength tensors and their covariant derivatives are the same as in Ref. [27]. In addition to the 16 CP-even operators, we also encounter 5 CP-odd operators, which we define as follows: with B µν , W a µν , and G A µν defined as in Eq. (43). We now discuss the SU(2) gauge traces appearing in our computation. The building blocks of these traces are the Higgs doublet H, its conjugate H = iσ 2 H * , and their covariant derivatives, which read: where the hypercharge of the Higgs doublet (conjugate Higgs doublet) is

C Relation to the UOLEA Coefficients
We now provide the relations between our coefficients, g N , to the (symmetrized) bosonic UOLEA coefficientsf N reported in Ref. [13]. As explained at the end of Sec. 2, the universal coefficients corresponding to operators containing P insertions are not captured by the bosonic UOLEA. Consequently, g i 10 , g ij 18 and g ij 19 will not appear below. To simplify the expressions, we define the sum of two masses as: Using this definition, we have: