Convergent Bayesian Global Fits of 4D Composite Higgs Models

Models in which the Higgs boson is a composite pseudo-Nambu-Goldstone boson offer attractive solutions to the Higgs mass naturalness problem. We consider three such models based on the minimal $SO(5) \rightarrow SO(4)$ symmetry breaking pattern, and perform convergent global fits on the models under a Bayesian framework in order to find the regions of their parameter spaces that best fit a wide range of constraints, including recent Higgs measurements. We use a novel technique to analyse the fine-tuning of the models, quantifying the tuning as the Kullback-Leibler divergence from the prior to the posterior probability on the parameter space. Each model is found to be able to satisfy all constraints at the $3\sigma$ level simultaneously. As a by-product of the fits, we analyse the collider phenomenology of our models in these viable regions. In two of the three models, we find that the $g g \rightarrow H \rightarrow \gamma \gamma$ cross section is less than ${\sim}90$% that predicted by the SM, which is already in slight tension with experiment and could potentially be ruled out in the future high-luminosity run of the LHC. In addition, the lightest fermions $F$ arising from the new strong dynamics in these models are seen in general to lie above ${\sim}1.1$ TeV, with the $F \rightarrow tW^{+}$ and $F \rightarrow \bar{b}W^{+}$ decays offering particularly promising channels for probing these models in future collider searches.


Introduction
The Standard Model (SM), assuming UV completion at a high energy scale, is only able to accommodate the observed mass of the Higgs boson through an excessive fine-tuning of the model parameters on account of the large quantum loop contributions the mass receives. Composite Higgs Models (CHMs), on the other hand, contend that the Higgs boson is in fact a bound state of some new strong dynamics whose confinement scale cuts off the loop contributions at O(1) TeV, alleviating this issue of fine-tuning [1][2][3]. In such models, the Higgs emerges as a naturally light pseudo-Nambu-Goldstone boson (pNGB) of some spontaneously broken symmetry of the new "composite" sector to explain the absence of other composite resonances in collider searches to date. Of considerable interest is the Minimal CHM (MCHM): the simplest viable CHM, based on the symmetry breaking pattern SO(5) → SO(4) that delivers only the Higgs doublet in the NGB spectrum while also ensuring custodial symmetry [4,5], which we focus on here. Many particular realisations of the MCHM exist, differing in their composite field content and their precise symmetry structures. It has proven particularly fruitful to consider those MCHMs that are holographically dual to certain 5D Gauge-Higgs-Unification theories, for they inherit calculability of the Higgs potential from the higher-dimensional models, and therefore possess extra predictivity in comparison to generic MCHMs [4][5][6][7][8][9][10][11]. A dimensional deconstruction [12,13] of the extra dimension leads to "elementary" boundary fields and a tower of "composite" resonances associated with the discrete points inside the extra dimension. The elementary fields mirror the non-Higgs SM fields and mix with the composite fields so that physical particles are "partially composite" superpositions of both [14]. We will be working with low-energy effective descriptions of such models known as Minimal 4D Composite Higgs Models (M4DCHMs), wherein the composite resonance towers are truncated to only a finite number of levels but calculability of the Higgs potential is still maintained [15][16][17]. See also Ref. [18] for a review of 4DCHMs.
Models of this form, which may still differ in the SO(5) representations of the composite fields, are popular choices in studies of composite Higgs scenarios on account of their predictivity and relative simplicity 1 . Early studies mostly examined the theoretical form 1 For studies of 4DCHMs with non-minimal symmetry structures, see for example Refs. [19][20][21][22].
-1 -of the Higgs potential and its consequences on the Higgs mass and fine-tuning, at first only including composite partners for the top quark in the fundamental representation [15], and later also including bottom quark partners and other representations, with some simple numerical scans performed to verify the results [23]. Composite lepton partners have also been considered, and have been found to give interesting fine-tuning effects [24,25]. More comprehensive scans of M4DCHMs have since been performed, for example: • to study the Higgs phenomenology in regions of parameter space that reproduce the correct electroweak (EW) scale and accurate quark and Higgs masses, for a wide variety of models with top and bottom quark partners in different representations [26], • to analyse points consistent with EW precision tests, Higgs physics, flavour physics, and (Run 1) LHC resonance bounds in four M4DCHMs that include all quark partners in the fundamental representation, subject to different flavour symmetries [27], • and to investigate the fine-tuning behaviour in regions fitted to give accurate EW scales and SM masses of models with composite partners for the top quark [28], and also for all third generation fermions [25].
However, no work so far has contained rigorous, convergent statistical global fits of such models. This is primarily on account of their large parameter spaces, making exhaustive parameter scans prohibitively computationally expensive, along with the fact that the parameters all have highly non-linear effects on both the SM and composite sectors, so that even simplified analyses of the parameter spaces are difficult. The purpose of the present work is to extend this research by providing the first full convergent global fits of several realistic M4DCHMs, mapping the regions of their parameter spaces that best fit a wide range of constraints. Given the difficulties mentioned above, we consider limiting models in which out of the SM fermions, only the top and bottom quarks couple to the composite sector, so that we may have more manageable parameter spaces. We will explore three different M4DCHMs, distinguished by the representations of the quark partners under SO (5). Our fits are performed under a Bayesian framework and are facilitated by a sophisticated nested sampling algorithm that delivers convergent results where other techniques have struggled. This Bayesian approach is particularly relevant to the case of CHMs, or indeed any proposed solution to the Higgs mass naturalness problem, for it takes into account both the experimental fitness and the naturalness of the models. We use the numerical techniques of Refs. [27,29] to calculate the observables at each point, and employ their same constraints (with updated values) where applicable. We also analyse the composite resonance cross sections and Higgs signal strengths predicted by the models in their viable regions to help guide future collider searches for new physics.
The structure of this paper is as follows. In Section 2 we give a brief overview of the M4DCHM and specify all three of our chosen fermion sectors, and in Section 3 we detail our procedure for scanning these models. We present and discuss our fit results, as well as make some comparisons between the different models based on their fits, in Section 4. We analyse experimental signatures of the three models in Section 5, and conclude with a summary of our findings in Section 6.

Model description
Here we outline the structure of the M4DCHM, first described in Ref. [15]. Throughout this work we refer specifically to the two-site M4DCHM, which contains only one level of composite resonances, for this is the minimum matter content needed for calculability of the Higgs potential. Site 0 consists of massless elementary fields with the same quantum numbers as the non-Higgs SM fields, while Site 1 contains the composite fields, which mix with the elementary fields both linearly and through proto-Yukawa terms to enact partial compositeness. An SO(5) → SO(4) symmetry breaking gives rise to the Higgs doublet, though the SO(5) symmetry is explicitly broken by the elementary-composite mixing so that the Higgs is really a pNGB. This mixing necessitates SU (3) and U (1) X symmetries on each site so that particles can be assigned the correct colour charge and the correct hypercharge under the definition where the isospin operator T 3 R is the third generator 2 of SU (2) R ⊂ SO (5). The overall symmetry breaking pattern is then (2. 2) The symmetry groups G on each site spontaneously break to their diagonal subgroup, leading to non-linear σ-model "link" fields Ω i (i = 1, 2) that parameterise the resulting NGBs. Most of them are would-be NGBs that become the longitudinal polarisations of the gauge fields on each site such that in addition to the SM bosons, there are ten massive SO (5) vector bosons, eight massive SU (3) resonances, and a massive abelian U (1) X resonance. The remaining physical (p)NGBs, parameterised by the product Ω = Ω 1 Ω 2 , constitute the composite Higgs field.

Boson sector
The boson sector of the theory is quite simple, for it is composed entirely of gauge bosons and NGBs. The low-energy dynamics of NGBs are entirely specified by their coset structure and are governed by the σ-model Lagrangian, so that the complete bosonic Lagrangian can be expressed as 3) The gauge sector consists of G 0 µ , W 0 µ , and B 0 µ -the elementary counterparts of the SM gauge fields -and ρ Gµ , ρ µ , and ρ Xµ -the SU (3), SO (5), and U (1) X gauge fields of the composite sector. They have the usual kinetic terms . (2.5) In order to give the σ-model Lagrangian, it will be useful to recast the link fields Ω i into independent link fields for each simple group factor in the global symmetry breakings. In this set is Ω 1,X,G respectively associated to the SO(5), U (1) X , and SU (3) c diagonal subgroup breakings between the two sites, and also Ω 2 for the breaking SO(5) → SO(4). Each will carry its own NGB decay constant f 1,2,X,G that determines the scale of the associated symmetry breaking. The link fields transform under the symmetries at each site as SO(5) 0 × SO(5) 1 : (4) : where transformations g k come from Site k, and h ∈ SO(4). Accordingly, their covariant derivatives are given by Here there are elementary gauge couplings g 0 , g 0 , and g 0 s mirroring the SM couplings, as well as new composite couplings g ρ , g X , and g G of SO(5) 1 , U (1) 1 X , and SU (3) 1 c . The generators T A = {T a L , T a R ,T a } of SO(5) have been grouped into SU (2) L , SU (2) R , and coset components (see Appendix A), allowing the SO(5) resonances to be expressed as Note that since SO(5) is spontaneously broken to SO(4), there is a vacuum vector Φ 0 left invariant under SO(4), so the product Φ := Ω 2 Φ 0 transforms as a fundamental of SO(5) 1 . We will use a basis where Φ 0 = (0, 0, 0, 0, 1) . With these, the leading-order σ-model Lagrangian is We now go to the Site 0 holographic gauge to remove the unphysical NGBs, setting Now the gauge boson mixing terms in Equation (2.9) become readily apparent: Hidden inside the first line of Equation (2.12) is mixing between a 4 µ and the Higgs: This mixing is removed, keeping the Higgs kinetic term conventionally normalised, with the field redefinitions (2.14) Note that the Higgs interactions now depend only on the quantity Hence, masses can only depend on the Higgs field through the "vacuum misalignment" parameter s h . We will return to this point in Section 2.5.

Fermion sector
At this point, the remaining freedom in specifying a particular model lies in establishing the fermionic content of the composite site and the representations of SO(5)×U (1) X under which it transforms. We will be restricting to the case where the composite sector mixes only with the third generation elementary quarks in the interests of reducing the parameter space of the theory. This is a reasonable limiting case, for the lighter fermions must mix with the composite sector only weakly if their left-and right-handed chiralities mix with roughly equal strengths. We will be considering three different models in this work: the M4DCHM 5−5−5 , the M4DCHM 14−14−10 , and the M4DCHM 14−1−10 , all of which have been outlined in Ref. [26]. The labels M4DCHM q−t−b specify the representations (q, t, b) of SO (5) in which the respective composite partners (Ψ q ,Ψ t ,Ψ b ) of the elementary (q 0 The structure of the elementary-composite interactions is depicted graphically in Figure 1. L actually couples to two different composite multiplets Ψ t and Ψ b , which respectively couple toΨ t andΨ b with independent interaction strengths. These models are chosen primarily because they provide custodial protection for the Zb LbL coupling from tree-level corrections [4], and because of their different relationships with fine-tuning. The fine-tuning of a model is commonly estimated as f 2 /v 2 , where v ≈ 246 GeV is the EWSB scale, though the M4DCHM 5−5−5 is subject to a parametrically higher "double tuning" because it requires precise cancellations of leading terms in the Higgs potential to achieve a non-trivial minimum and to produce the correct EWSB scale [23,30]. The other models need not feature such extra tunings on account of extra leading-order invariants possible in the Higgs potential with the 14 representation. The M4DCHM 14−14−10 and M4DCHM 14−1−10 each also offer interesting prospects, with the former having symmetries that allow for two independent proto-Yukawa couplings for the top quark, and the latter being able to feature an entirely composite right-handed top quark (though we do not take it as composite in this work). These are not the only possible models that fit our criteria, but it has been found in Ref. [26] that the various other representation combinations either tend to give experimentally unfavourable predictions, or are qualitatively similar to the three we consider here.
In specifying the Lagrangian for the fermion sector of each model, we split up the contributions as The high-energy Lagrangian for the (partially) composite quarks, consisting of the third generation elementary quarks and their composite partners Ψ andΨ, will be given below for each model. Elementary fields will be embedded into incomplete multiplets ψ in the same representation as their partners for convenience in writing gauge-invariant interactions 3 .
3 A more general approach is to have the SM-like elementary fields embedded in complete representations, with mass terms given to the supplementary fields in the elementary multiplets to explicitly break the SO(5) symmetry [31,32]. Our formulation is obtained in the limit where these mass couplings go to infinity. Along with the assumption of maximal symmetry [33], such general models alleviate issues with light quark partners and promise favourable levels of fine-tuning. It would be interesting to fit these models in future work.
-6 -In all cases the covariant derivatives of the fields are given by , and λ a are the SU (3) generators, taken so that each elementary or composite field couples with the same strength to the SU (3) gauge bosons. It should be kept in mind that a multiplet Ψ in the symmetric 14 or antisymmetric 10 of SO(5), when expressed as a matrix, is acted upon by the generators as (2.18) Refer to Appendix A for the explicit embeddings of fields into representations of SO(5).

M4DCHM 5−5−5
Unlike the other models we will be exploring, the M4DCHM 5−5−5 requires two composite multiplets Ψ t,b to couple to q 0 L , simply because a single partner multiplet does not provide sufficient couplings for both the top and the bottom quark to have mass. The partner content of this model consists of Ψ t ,Ψ t in the 5 + 2 3 representation of SO(5) × U (1) X , and . Their interactions are specified by the quark Lagrangian Here the third generation elementary quarks are embedded into the incomplete SO(5) fundamentals Note that this is not the most general Lagrangian obeying the symmetries of the theory. Terms such as ψ L Ω 1ΨR orΨ L Ψ R , if included, would reduce the model to the two-site Discrete CHM and lead to a logarithmic divergence in the Higgs potential. Similar considerations exist for the other representations.

M4DCHM 14−14−10
Here we have two composite multiplets, Ψ q andΨ t , in the traceless symmetric 14 of SO(5), and one multipletΨ b in the antisymmetric 10. Each composite multiplet carries an X charge of + 2 3 . All of them transform adjointly, so invariant terms may be found by taking the traces of various field combinations. Writing out the invariant interactions yields the quark Lagrangian The elementary quark fields are embedded into incomplete representations as Notice that by the symmetries of this model there is an extra allowed Yukawa-like term, having a coupling strengthỸ t . The analogous terms proportional to m Y b andỸ b vanish due to the symmetry (antisymmetry) of Ψ q (Ψ b ).

M4DCHM 14−1−10
The M4DCHM 14−1−10 has a similar structure to the M4DCHM 14−14−10 , except withΨ t now being a single field in the 1 + 2 3 of SO(5) × U (1) X instead of a multiplet. The quark Lagrangian is given by The elementary ψ q L and ψ b R multiplets are the same here as in Equation (2.23), but now ψ t R = t 0 R . Notice that the symmetries do not allow any off-diagonal mass terms m Y in this model.

Elementary fermions
In a truly consistent CHM, all of the lighter fermions would be incorporated through partial compositeness in a similar manner to the above. However, in our work, we treat them as entirely elementary for simplicity. We reproduce the SM lepton sector with bilinear couplings to the low-energy Higgs pNGB: with no Higgs couplings, and again coupling to the elementary gauge fields. Admittedly there is no reason to treat the mass terms of the elementary quarks and leptons differently. However, the difference is largely inconsequential due to the small couplings to the Higgs involved, and because there is little distinction between a mass term and a Higgs coupling below the EW scale.

Higgs potential
Mixing between the elementary and composite sites explicitly breaks the global symmetries of the composite sector, and therefore leads to a quantum effective potential for the (now pseudo-) NGB Higgs. At the low energies that we are capable of probing, it is useful to integrate out the composite degrees of freedom and work with the effective (composite quark) Lagrangians that, in momentum space, will be of the form for some model-dependent correlators Π ψ and M ψ provided in Appendix B. With these, the fermion contributions to the effective potential, with only partially composite third generation quarks, can be shown to be where N c = 3 is the number of colours and a Wick rotation to Euclidean space was used. Notice that scaling the mass parameters does not change the stationary points of the potential. The location of its minimum -the Higgs VEV h -can be found by expanding the potential in powers of s h as V eff (h) =: −γs 2 h + βs 4 h to obtain the coefficients γ and β that dictate the VEV satisfies and which can further be used to calculate the Higgs mass as We neglect detailing the incorporation of the gauge boson contributions to the potential in this framework, since they are of secondary importance. Suffice to say that the gauge boson contribution to γ can be calculated analytically, with the result that [25] to first order in t θ := g 0 /g ρ . Here we have the approximate masses of the lightest composite gauge bosons (see Section 2.5) We use this formula in our scans, and neglect the contributions to β from the gauge bosons.

Particle content
One interesting signature of composite Higgs models is the existence of new resonances in the few-TeV mass range. In this section we give an overview of the particle content in each of the models we are considering, including some analysis of the particles' masses.

Boson sector
Each model shares the same boson sector. It contains, in addition to the SM bosons, the heavy gluons, five other neutral bosons, and three bosons of unit electric charge. Their mass mixing matrices are given as functions of the Lagrangian parameters in Appendix C. Singular values of the matrices are the squared (tree-level) masses of the resonances, though not all have useful analytic forms. The bosons with easily calculable masses are the heavy gluons and the neutral a 4 resonance, which have squared masses of 33) and the charged bosons, whose mass spectrum is given by m 2 charged = 1 2 f 2 1 g 2 ρ , m 2 1 , m 2 2 , m 2 3 , where the last three masses are the solutions to The W boson mass, for example, distinguished as being O(s h ), is given by Neglecting the higher-order terms, we can match to the SM expression m W = 1 2 gv to find the relations Similar matchings can also be done to yield relations for the remaining SM gauge couplings in terms of the elementary and composite couplings that we use throughout this work: , and 1 g 2 (2.37) The same sort of analysis can be done for the other resonances, but the contributions to their masses from the Higgs field are only relatively minor. We simply give their masses in the s h → 0 limit, which will be approximately satisfied by realistic points and can be used as a zeroth order estimate. In this limit, excluding the electroweak bosons (which will be massless), the charged bosons have masses Here we have defined the functions which will also be useful in expressing the fermion masses. Note there is some degeneracy among the masses of the charged and neutral bosons in the s h → 0 limit. In fact, numerical calculation of the masses for realistic parameter points with small non-zero s h suggests that two charged bosons will always separately be very close in mass to two neutral bosons, and the last charged boson will be very close in mass to two other neutral bosons ("very close" here meaning within ∼0.1%).

Fermion sector
The fermion sector of each model can be analysed using the explicit embeddings of fields in the different representations of SO(5) given in Appendix A. We group the particles according to their electric charge, which is determined through the relation   Each icon corresponds to one particle, or to a pair of particles for masses given in terms of the functions M ± defined by Equation (2.40). Some masses cannot be expressed in a reasonably understandable form, in which case they are left unspecified.
From the fields' SU (2) L × SU (2) R quantum numbers and the U (1) X charges of the multiplets specified in Section 2.2, the models are seen to deliver some up-type (U) and down-type (D) resonances, along with particles of exotic charge Q 4/3 , Q 5/3 , and Q 8/3 . The number of new particles in each model is listed in Table 1.
The (tree-level) masses of the fermions are found as the singular values of the fermion mass matrices in Appendix C. For convenience, the non-SM fermion masses in each model are provided in Figure 2. Similarly to the bosons, many fermion masses are not easily expressible analytically, and we resort to giving some masses in the s h → 0 limit. Note the many degeneracies and approximate degeneracies among the particle masses in each model.
- 12 -In this section we provide details regarding our global fits of the models specified in Section 2. A description of the scanning algorithm used for the fits is given in Section 3.1, and our treatment of the scan parameters in Section 3.2. The experimental constraints employed to determine the validity of a given parameter point are outlined in Section 3.3.

Scanning algorithm
The purpose of our fits is twofold: firstly, to obtain estimates of the parameters of each model given current experimental data, and secondly to facilitate model comparisons under a Bayesian framework through the calculation of their Bayesian evidences. Both are made possible by scanning the parameter spaces with the nested sampling algorithm [34], which has proven effective at sampling difficult spaces in many physical theories, including CHMs. For our fits we make use of the recently developed nested sampling package PolyChord, which uses a sophisticated version of slice sampling to generate points efficiently [35,36].
Nested sampling algorithms explore a model's parameter space using a certain number n live of "live" parameter points p that are iteratively constrained into progressively smaller regions in which the predictions of the model better match observed data. The experimental constraints are enacted through a likelihood function L(p) that the algorithm tries to minimise. We take the likelihood to be Gaussian in the observables: where O theo (p) is a "vector" of the predicted values of the observables that have experimental values O exp , and C is the covariance matrix that takes into account theoretical and experimental uncertainties, and also correlations between observables. Details about the specific data we use are given in Section 3.3. The live points are sampled according to a prior distribution π imposed on the space, which together with the likelihood L gives the posterior probability P of each point is the Bayesian evidence of the model in question. Throughout nested sampling exploration, the Bayesian evidence is able to be estimated from the prior-weighted volumes of the live points at each iteration and their likelihoods. We declare the algorithm "converged" when the evidence carried by the remaining live points falls below 10 −3 times that contributed from the previous iterations. The posterior distribution and the evidence both weigh the fitness of a point or region of space against its naturalness, and so are well suited for judging the viability of composite Dimensionality 19 17 15 Table 2: Parameters present in each model.
Higgs models. A similarly helpful quantity that may be defined is the Kullback-Leibler (KL) divergence between the posterior and the prior, which measures, in a sense, the extra information gained when going from the prior to the posterior [37], and therefore is an indirect measure of the fine-tuning of a model. Using Equation (3.2), this may be rewritten in terms of the average loglikelihood over the posterior as [38,39] ln(L) P = D KL + ln(Z). (3.5) Once the evidence and KL divergence are known, this is a useful way to calculate the posterior-averaged log-likelihood as a quick check of how well the model fits the data.

Scan parameters
A full list of the parameters for each of our models is given in Table 2. We take the approach of Refs. [25,26,28] and do not scan over these parameters as presented, but instead scan over the mass-dimension parameters normalised by the scale f . This choice is convenient because the Higgs potential can be minimised 4 to find the misalignment s h , and the mass scale only set later by defining in accordance with Equation (2.36) to automatically reproduce the correct EWSB scale, and by extension the experimental W and Z boson masses. Our treatment of the parameters in our scans is summarised in Table 3, and explained in greater detail below.

NGB decay constants
The Higgs decay constant f is one of the more consequential parameters in the theory, defining the approximate mass scale m * ∼ g ρ f as well as the naïve energy cutoff Λ f = 4πf . The ranges of other parameters are limited by the value of f . For example, f 1 must be greater than f by virtue of Equation (2.14), and less than √ 3f to maintain partial unitarisation of NGB scattering [17]. The other decay constants are constrained by to avoid decoupling any resonances. We use these constraints to define the bounds for these parameters, with the added condition that all are greater than 0.5 TeV. Points that are generated within the bounds but that do not satisfy these consistency conditions are immediately discarded -a process that simply amounts to a modification of the prior such that only consistent points are admitted.

Gauge couplings
Given that the composite sector is strongly coupled and we can only perform calculations in the semi-perturbative regime, we must take the gauge couplings in Table 2 to be between 1 and 4π. In addition to these couplings, we also vary the SM gauge couplings within their experimental limits (at the scale of the top mass). The elementary gauge couplings are then determined by the relations in Equations (2.36) and (2.37). For g 0 s to be real we require that g G > g s . We also impose the restrictions to avoid vector resonance masses above the cutoff Λ f .

Link, mass, and proto-Yukawa couplings
The remaining parameters in Table 2 have dimensions of mass and so are given upper bounds of Λ f . Instead of directly scanning over the off-diagonal masses and proto-Yukawa couplings in Table 2, however, it turns out to be more convenient to scan over particular linear combinations that appear in the correlators and mass matrices. Specifically, we scan over normalised by f . All of these parameters that we scan over are taken to be positive through field redefinitions, and in our initial test scans they were given arbitrarily chosen lower bounds of e −8.5 . These ranges were further restricted into those listed in Table 3 based on the results of these tests. It was also apparent that the experimental constraints were drawing many of these parameters towards lower values, so scanning them with logarithmic priors was found to be advantageous.  Table 3: Ranges and priors used for the parameters in our scans. The masses m ρ and m a , defined in Equation (2.32), reparameterise the decay constants into a form more suitable for calculating the gauge boson contribution to the Higgs potential. The normalisation factor f is determined after the potential is minimised. Some parameters are subject to further constraints specified in the main text.

Constraints
We employ a wide range of observables to constrain the models, consisting of SM masses, EW radiative corrections, Z boson decay ratios, and Higgs signal strengths. Bounds on the production of new heavy resonances from direct collider searches are also considered for realistic points found in the scans, though we chose not to include these constraints in the scans because of the increase in computing time they would require. We use software developed for Ref. [27] (and the subsequent modifications for Ref. [29]) to calculate the observables and likelihood for each point, so we only give a brief overview of each constraint below and refer the reader to Section 3.1 of Ref. [27] for details on their implementation and further discussion. Only those constraints applicable to our third-generation-partiallycomposite-only models are included. The experimental values have been updated for this work and are given in Table 4.
m t (173.0(4) ± 1.0 theory ) GeV [40,41] Table 4: Experimental values used for constraints in our scans. The observables are grouped (in order) into SM masses, oblique parameters, Z decay ratios, and Higgs signal strengths. The latter employ multiple independent measurements. Note the values for S and T have a correlation coefficient of +0.92 that is taken into account in the χ 2 calculation.

SM masses
All of the masses of SM particles that are predictions of the theory -namely, the top and bottom quark masses and the Higgs mass -are used as constraints. Once the Higgs mass is found by the method outlined in Section 2.4, the other masses are found at tree level by -17 -diagonalising the mass matrices in Appendix C at the point h = h . The top and bottom quarks are identified as the third lightest up-and down-type particles in the theory.

Oblique parameters
Important constraints for CHMs come from EW precision observables, which restrict the non-linear dynamics of the pNGB Higgs and place lower bounds on composite vector resonance masses. Such observables are conveniently parameterised by the Peskin-Takeuchi S and T parameters [47,48], which we take as constraints. The S parameter is calculated at tree level and mostly constrains the vector resonance masses. The T parameter arises at one-loop and is given in terms of the vacuum polarisation of the electroweak bosons. We consider only the dominant contributions from fermion loops, including in particular the effects of composite fermion resonances. Absolute theoretical uncertainties of 0.05 and 0.10 are assigned to the respective theory predictions of S and T and they are assumed to be uncorrelated in contrast to the experimental uncertainties.

Z decays
The compositeness of the third generation quarks modifies the Z boson couplings, and in particular Zb LbL . These couplings are tightly constrained by the Z boson decay widths and we consider their ratios is the total Z hadronic width.

Higgs signal strengths
We use Higgs signal strengths from gluon fusion production to further constrain the nonlinear dynamics of the pNGB Higgs. The signal strength µ gg X for the decay into a final state X is defined as the ratio of the measured cross section to the predicted SM value: We include as observables µ gg X for X = τ τ , W W , ZZ, and γγ. The partial widths for decays into massive particles are calculated at tree level, and those for decays into massless particles are calculated at one-loop order including all fermionic and bosonic contributions. For each of these observables we use three independent measurements: one using combined ATLAS and CMS data from Run 1 of the LHC, and the others separately from ATLAS and CMS using Run 2 data. The effects of these constraints on the parameter spaces of M4DCHMs and similar models have recently been analysed in detail in Refs. [49,50].

Collider searches
Searches for heavy resonances in colliders place useful upper bounds on the production cross section times branching ratio of composite resonances for various decay modes. As mentioned above, we do not include these as constraints in the scans themselves, and only analyse them for realistic points found in the scans after convergence. All of the experimental analyses we use are featured in Tables 7 to 11, in Appendix D. This builds on the list of analyses used in Ref. [29] with an additional 40 LHC searches at √ s = 13 TeV.
Note that since the searches do not provide measured values and uncertainties but upper bounds, there is a slight departure from Equation (3.1) in calculating the χ 2 contribution of each decay. In this case, the 95% CL upper bound on the cross section times the branching ratio, denoted by X bound , is taken into account in such a way that a theory prediction X theo of the same value should correspond to a χ 2 = 1.96 2 , i.e. the χ 2 contribution of a particle with X theo in a given decay channel is where σ X = X bound /1.96. While this is applicable for the case that the expected bound X exp equals the observed bound X obs , we correct this relation in the case that both are different by using σ X = X exp /1.96 and defining The function κ corrects the χ 2 contribution if X obs /X exp = 1 while for X obs /X exp = 1, we have κ(1) = 1, i.e. in this case we have x = 0 and Equation (3.14) reduces to Equation (3.13).
In calculating the χ 2 for each decay, we only use the maximum χ 2 contribution from all the analyses restricting the decay. While a tighter bound would be given by combining the contributions from independent analyses, we have found that the improvement in the excluded cross section at the 3σ level when doing this is often smaller than the error introduced by the narrow width approximation used in our calculations, and so implementing this would not be worthwhile.

Mass cuts
As a final measure, we impose harsh penalties for any new fermionic resonances below 500 GeV because such a resonance would likely have already been discovered. Ruling out points using hard mass cuts below 500 GeV makes it prohibitively time consuming to find an initial population of viable points, so instead we assign a steep one-sided Gaussian likelihood to the fermionic masses. This makes initial populations of points easier to find, and they will evolve towards points that give more reasonable masses as the scans proceed.

Results
The results of our global fits of the M4DCHM 5−5−5 , the M4DCHM 14−14−10 , and the M4DCHM 14−1−10 are presented below in Sections 4.1 to 4.3, and are collectively discussed -19 -in Section 4.4. As a by-product of the fits, we analyse the experimental signatures of the three models in Section 5 in the hope of guiding the search for evidence of these models at the LHC.
For each of the models, we performed multiple scans to verify the reproducibility of the results. This was a lengthy process, with each of our production-ready scans taking at least one week to converge running on 64 cores, or equivalent. The levels of agreement between the scans are shown in Appendix E. Broadly, results for the M4DCHM 14−14−10 and the M4DCHM 14−1−10 are highly reproducible, while there is somewhat less agreement in the results for the M4DCHM 5−5−5 . We analyse the combined samples from all scans for each model 5 using anesthetic [51] in order to maximise the robustness of our conclusions. In the following discussion we will focus on only those features of the results that are seen across multiple scans. We look first at the results for the gauge parameters, in Figure 3. There is a clear pattern for the NGB decay constants: their priors all favour lower values of the constants, successfully encoding the notion that smaller decay constants are more natural (cf. discussion in Section 2.2), while their posteriors occupy the higher values, indicating the smaller decay constants are experimentally disfavoured as is already known from precision EW tests [7]. The posterior for f , for example, lies mostly between 1.2 TeV and 1.9 TeV, with a maximum at ∼1.4 TeV. The upper bound here comes from suppression from the prior, while the lower bound is due to a variety of constraints; namely, the oblique constraints, the Higgs signal strength constraints, and the Z decay constraints. It is not impossible for such constraints to be individually well-satisfied at low f ; it is simply that many points in this range do not satisfy the constraints well, and this large low-likelihood region suppresses the posterior over this range because of the marginalisation process 6 .
There does not seem to be any structure in the posteriors of the other decay constants in Figure 3 beyond the consistency conditions outlined in Section 3.2. That is, the SO(5) decay constant f 1 has no particular preference for any values between its theoretical bounds of f and √ 3f , having a posterior peaking around 2 TeV, and likewise for the U (1) X and SU (3) c decay constants f X and f G , which have similar posterior distributions between ∼1 TeV and ∼4 TeV, though the results for f G are not particularly consistent. These parameters are constrained by the same factors as f mentioned above, and this is true across all models.
As for the gauge couplings, the posterior for the SO(5) 1 coupling g ρ extends across the entire range [1, 4π] with a preference for values between ∼2 and ∼8, driven largely by the prior with some additional experimental factors at play, while no definitive results were found for the U (1) 1 X and SU (3) 1 3 couplings g X and g G . The Higgs signal strength, SM mass, and oblique constraints all seem to slightly disfavour values for g ρ at the extreme ends of the range, while the Z decay constraint R b disfavours values below ∼2, in opposition to R e and R µ , which favour those values. It is unsurprising that there is little structure in the -21 - posteriors of f G and g G since we do not include any constraints on heavy gluon decays or flavour physics in our scans. For the top sector parameters in Figure 4, the posteriors are significantly localised in comparison to the priors, demonstrating the effects of the constraints on the space. However, the posteriors all peak around the same regions as the priors, so it seems that satisfying the constraints does not require excessive fine-tuning. Indeed, the elementarycomposite mixing and composite mass parameters all tend to be constrained around ∆ tL ∼ ∆ tR ∼ m t ∼ mt ∼ 1.8 TeV, which, from the findings of Ref. [23], indicates that double tuning is not so much of an issue. All of these parameters are primarily constrained by the SM masses, except for the lower bound of 0.5 TeV on mt that comes from the resonance -22 - mass cutoff 7 . Not particularly well constrained is the off-diagonal mass m Yt , which may take any value below ∼2 TeV. The last top sector parameter, m Yt + Y t , clearly prefers larger values from around 2 TeV to 20 TeV, and this is due largely to the SM mass and oblique constraints. Finally, we discuss the bottom sector parameters in Figure 5. There is not as much uniformity among these as for the top parameters, presumably because these play a subdominant role in the Higgs potential. Interestingly, ∆ bL is very narrowly constrained, between around 10 GeV and 30 GeV. Much of this is due to the prior, which suppresses values above ∼150 GeV, but again there are other factors present. The lower bound on its posterior is easily identified as coming from the bottom quark mass constraint, while the upper bound comes more subtly from a large region of points with ∆ bL 50 GeV that significantly violate the Z decay constraints. The right-handed bottom quark coupling ∆ bR , on the other hand, is quite large, ranging between roughly 1 TeV and 10 TeV, with the lower bound coming from the bottom quark and Higgs mass constraints. The bottom mass constrains mb to lie above ∼1 TeV, and the prior below ∼5 TeV. The final parameter, m Y b + Y b , lies mostly below 2 TeV, with a preference towards the lower values stemming from the SM mass and Z decay constraints. Results for m b and m Y b are not particularly consistent, but both generally occupy large values, with the former ranging up to ∼4.5 TeV, and the latter up to ∼10 TeV.
One interesting feature to note here is that the posteriors of m b and m Y b display a strong correlation that is not present in the prior, and therefore must be a result of one or more of the experimental constraints, though it is difficult to pinpoint the exact cause. Other structures to note in Figures 4 and 5 are slight correlations in the priors and posteriors of the pairs (m t , m Yt + Y t ) and (m t , mt), a mild correlation in the posteriors of (∆ bR , mb), and posteriors that seem to avoid small values for both parameters in the pairs (∆ tL , m t ), (∆ tL , m Yt + Y t ), and (∆ bR , m Y b ). More will be said about these results when comparing models in Section 4.4.

M4DCHM 14−14−10
Results for this model are given in Figures 6 and 7. The priors and posteriors for the gauge sector parameters, in Figure 6, look much the same as in the M4DCHM 5−5−5 , with the NGB decay constant posteriors concentrated well above the small values preferred by the prior. The posterior for f enjoys a slightly higher range than before, from around 1.5 TeV to 2.2 TeV. Tending towards the higher end of its range is f 1 , spanning between ∼1.5 TeV and ∼3.8 TeV, while the remaining decay constants are constrained as 1.2 TeV f X 7.5 TeV and 1.1 TeV f G 6.5 TeV. The gauge coupling g ρ has a similar distribution to before, preferring mid-range values, again seemingly driven by the prior in combination with the Higgs signal strength and Z decay constraints. As before, g X and g G were not sufficiently well constrained to give consistent results. However, it seems that stronger g G couplings are disfavoured as they appear to entail undesirably small values of f G .
The M4DCHM 14−14−10 fermion parameters are seen from Figure 7 to be significantly more constrained than those in the M4DCHM 5−5−5 , with posteriors concentrated on considerably smaller regions in comparison to the priors -a feature emblematic of higher fine-tuning. This is particularly evident for the on-diagonal mass parameter m t , whose prior allows it to take effectively any value up to 10 TeV but whose posterior is localised between ∼0.6 TeV and ∼2 TeV. Mass constraints are the main factor here, with the Z decay and Higgs signal strength constraints also contributing to the upper bound somewhat, and the same is true for all of the other on-diagonal mass and mixing parameters as well.
has a posterior constrained to lower values 0.8 TeV. In other words, it must be that Y t ≈ −2m Yt . It is not clear exactly which experiments are to blame, but it seems that -25 -   The gauge coupling g ρ behaves in the same way as in the other models, both in its preference for values from ∼2 to ∼8, and in its relationships to the various constraints. However, the same is not true for the other gauge couplings. Surprisingly, in this model the posteriors for g X and g G almost exactly align with the priors, and indeed there does not seem to be any significant dependence on these parameters for any of the observables.
In contrast to the gauge parameters, there is a clear difference in the posteriors of the fermion parameters between the M4DCHM 14−1−10 and the M4DCHM 14−14−10 : those in the M4DCHM 14−1−10 , shown in Figure 9, cover a markedly larger fraction of volume of the total parameter space, so these parameters require less fine-tuning to satisfy the constraints than those in the M4DCHM 14−14−10 . This is taken to the extreme by the mass parameter m t being almost entirely unconstrained in this model: it is able to take values up to ∼10 TeV but with great preference towards values 0.5 TeV purely stemming from the prior. Other parameters that require little fine-tuning, having posteriors well within the regions favoured by their priors, are ∆ q , ∆ b , and m q . The posteriors of the first two of these respectively span from ∼0.9 TeV to ∼5.0 TeV, and from ∼70 GeV to ∼1.8 TeV, again reflecting the right-handed bottom quark being mostly elementary. The last of these, m q , along with m b , has a hard lower bound of 0.5 TeV from the resonance mass cutoffs, with the former ranging upwards of 8 TeV and the latter being drawn up to ∼21 TeV by the constraints. If the LHC constraints were included in this fit, we might expect these lower bounds to be increased and larger parameter values to be favoured. Unlike the other models, here the right-handed top quark can mix very strongly with the composite sector, with the posterior for ∆ t going between ∼0.5 TeV and ∼17 TeV, being maximised at ∼7.5 TeV despite the prior giving preferable weight to lower values. Similarly, the proto-Yukawa couplings Y t and Y b are drawn to unnaturally large values by the constraints, with posteriors ranging from ∼1.6 TeV to ∼22 TeV, and from ∼50 GeV to ∼63 TeV respectively. The posteriors for all of these parameters are shaped primarily by the SM mass constraints, along with some contributions from the Z decay constraints that slightly disfavour very large values of the couplings. There is a slight positive correlation between the posteriors of m b and Y b , and a negative correlation in both the priors and posteriors of ∆ b and Y b . There is also a notable effect wherein large values of Y t are favoured in both the prior and posterior when Y b is also large. The prior tends to disfavour both ∆ q and Y t having small values, and to a lesser extent the same is true for both ∆ q and m q , and the posteriors for these pairs seem to be constrained -28 - to regions for which ∆ q Y t and ∆ q m q lie above some respective fixed values. Unfortunately, analysing the causes of these features would be too involved a process to be viable for us. As a final note, our result that ∆ t is usually much greater than m t indicates that the right-handed top quark must be almost completely composite in this model given the constraints. It might therefore do well to consider a model with the same 14 − 1 − 10 symmetry structure that incorporates t R as an entirely composite state as a limiting case of the M4DCHM 14−1−10 .

Discussion and model comparisons
Having seen the fit results for all of our models, we are now in a position to discuss the similarities and differences between the models, and compare them based on their fits. The Bayesian evidences Z calculated from the combined samples from the scans of each model are given in Table 5, along with other quantities of interest such as the KL divergence D KL , the posterior-averaged log-likelihood, and the maximum likelihoods found in the scans. The most obvious point to note from Table 5 is that the M4DCHM 5−5−5 has by far the greatest evidence of the three models, being ∼7800 times greater than the evidence for the M4DCHM 14−1−10 , itself ∼15% greater than that for the M4DCHM 14−14−10 . Granted we assign an equal prior likelihood to each of our models, the M4DCHM 5−5−5 would then be decisively preferred over the others, while there is no conclusive preference between the M4DCHM 14−1−10 and the M4DCHM 14−14−10 . But note that these evidences are priordependent; the choice of the prior distribution as well as the lower and upper bounds of the prior parameter ranges affect these results. A benefit of our sampling method is that it is possible to quantify the degree to which the evidences and other quantities in Table 5 are dependent on the prior distribution (but not parameter bounds) by re-weighting our samples as if a different prior distribution was chosen. Uniform priors on the fermion parameters, for example, lead to Table 6, where it is seen the evidences are considerably different in comparison with the logarithmic priors. While in this case the M4DCHM 5−5−5 still seems to be superior compared to the other models, this, however, also depends on the parameter bounds that have been determined in the initial test scans and which are listed in Table 3. In particular, for the left-handed composite-elementary mixings ∆ t L and ∆ b L in the M4DCHM 5−5−5 , we find that the bounds on ∆ t L are larger than those on ∆ b L by a -30 - factor of roughly the size of the top to bottom mass ratio. This means that in our priors for the M4DCHM 5−5−5 , the hierarchy among the third generation quark masses is present in the bounds of left-handed composite-elementary mixings, while this is not possible in the other two models, which have only a single left-handed mixing parameter ∆ q . The fact that different lower bounds have a much smaller effect in the case of uniform priors might be a reason why the evidences in Table 6 do not show such large differences between models as those in Table 5. However, also in this case, the prior parameter ranges are favourable for the M4DCHM 5−5−5 . Therefore, a fully conclusive model comparison would require a detailed analysis of the impact of prior distributions and bounds on the evidences and the other statistics shown in Table 5, and it would be interesting to perform such an analysis in future work.
Apart from the prior dependence, there are also other reasons for the hierarchy of preferences that can be ascertained from Table 5. First, the evidence for the M4DCHM 14−14−10 is relatively low primarily because of a large KL divergence, i.e. large fine-tuning. We have already recognised this tuning by way of the narrowly constrained posteriors in Figures 6  and 7. On the other hand, the M4DCHM 14−1−10 is disfavoured mostly because it has a much smaller posterior-averaged log-likelihood, i.e. because it has a more difficult time satisfying the experimental constraints. We find the Higgs signal strength constraints to be the main factor here. The M4DCHM 5−5−5 is the only model with both a low KL divergence (tuning) and relatively large average likelihood, which combine to give the greatest evidence. But here too, fine-tuning is wholly dependent on the prior distribution, as this defines the "natural" distribution of the parameters. The worsened tuning of the M4DCHM 5−5−5 under the uniform priors in Table 6 makes it clear that this model benefits greatly from the logarithmic priors on the fermion parameters in comparison to the other models. Note that fine-tuning has previously only been quantified in the literature of CHMs through the sensitivity of the observables with respect to the input parameters, either to first order with the Barbieri-Giudice (BG) measure [23,26,27,29], or with a higher-order -31 - measure [25,28,52], and our approach of using the KL divergence as a fine-tuning measure has not been taken before.
It is of interest to compare the above findings to the naïve fine-tuning estimates suggested by the distributions of the f parameter in each model. Recall the BG tuning required to reproduce the correct EWSB scale is ∼f 2 /v 2 , so the posterior distributions we have found for f correspond to relatively reasonable tunings of ∼25−60 in the M4DCHM 5−5−5 , ∼40−80 in the M4DCHM 14−14−10 , and ∼30−80 in the M4DCHM 14−1−10 . Note that despite its double tuning, the M4DCHM 5−5−5 has both the lowest KL divergence and the lowest naïve tuning of the three models, corroborating the findings of Ref. [26]. A further advantage of the M4DCHM 5−5−5 is apparent in the priors on f , displayed in Figure 10.
Here it is seen that f is generated at higher values more often in the M4DCHM 5−5−5 than in the other models, so the M4DCHM 5−5−5 can accommodate larger values of f with lower levels of tuning.
It also proves useful to analyse the priors for observables taken as constraints in each model, to see whether any model gives realistic predictions more naturally than others. We give the priors for the masses of the SM particles in Figure 11. As expected, the Higgs, being a pNGB, is generated most often with a low mass. But it also lies above 1 TeV relatively often, especially in the M4DCHM 14−14−10 , and to a lesser extent in the M4DCHM 14−1−10 , both of which have been noted in other works as having a tendency for producing too large a Higgs mass [23,26]. The top and bottom quark masses, on the other hand, all tend to be generated relatively close to their experimental values 8 . The M4DCHM 14−14−10 gives quite an unfavourable prior for the top mass, with very light values 80 GeV being preferentially generated. Remarkably, the prior for the top mass in the M4DCHM 5−5−5 has a maximum almost exactly at the experimental value. It should be stressed that this is not a feature of the model itself, but rather a consequence of the bounds imposed on the top -32 -sector parameters. Indeed, if the priors for the bottom sector and top sector parameters were switched in this model, then so too would the priors for the bottom quark and top quark masses switch. This is another example of the prior dependence discussed above which makes it easier for the M4DCHM 5−5−5 to reproduce the top mass, leading to a lower KL divergence and a larger evidence.
We conclude this section with some general comments on the parameter spaces of our models. During our scans, we noticed that the slice sampling steps take O(100) likelihood evaluations in each slice instead of the more typical O(5) evaluations, indicating that the likelihood is quite "sheety": it is only non-negligible across thin hypersurfaces in the parameter spaces. This notion can be quantified with the Bayesian model dimensionality (BMD) of each model, defined similarly to the KL divergence Equation (3.4), which measures the effective dimensionality of the region of parameter space that the posterior occupies [53]. This definition is used because it matches the actual dimensionality when the posterior is Gaussian. Broader posterior distributions will have lower BMDs, and narrower distributions higher BMDs. It can be difficult to intuit the BMD of a model based on the marginalised posteriors such as those given in Sections 4.1 to 4.3, since it is most often the case that the constraints restrict only certain non-trivial combinations of the input parameters, and these correlations would be obscured in a 2D plot by the marginalisation. In any case, these effective dimensionalities can be calculated with anesthetic and are listed for our models in Table 5. Clearly, the three models are indeed very sheety, with effective dimensionalities of ∼7 in the M4DCHM 5−5−5 and M4DCHM 14−1−10 , and ∼12 in the M4DCHM 14−14−10 , in contrast to their respective true dimensionalities of 22, 18, and 20 (including the SM gauge couplings). Positive constraints that restrict the space to reproduce particular values for the observables are evidently not sufficient to localise the parameters completely, and leave these ∼7 − 12 degrees of freedom within the spaces. These effective dimensionalities should be expected to be quite small on account of the fine-tuning required for EWSB and the effect of the constraints on the space, but it is difficult to say whether ∼7−12 is a reasonable number of dimensions on purely theoretical grounds.
It is further possible to calculate the principal directions in parameter space that point along these hypersurfaces with anesthetic. If these models were to be explored in future scans, it would be beneficial to use these directions to cleverly parameterise the spaces and scan along these optimal surfaces directly, to achieve greatly improved scanning efficiencies.

Experimental signatures
In our fits, we have found a number of points in each model that satisfy all constraints and direct collider search bounds at the 3σ level individually, which we regard as valid points. We have collected, in total, ∼14,000 valid points in the M4DCHM 5−5−5 , ∼77,000

Composite resonances
We begin with the BSM particle spectra. Recall that these models contain several uptype (U) and down-type (D) composite fermions, exotic fermions Q x with electric charges x ∈ {4/3, 5/3, 8/3}, as well as composite gauge bosons with unit charge (W) and neutral charge 9 (Z). Our convention for distinguishing the different members of each species is to label them with numerical subscripts in order from lightest to heaviest, including SM states. For example, the lightest neutral composite gauge boson will be called Z 3 , being the lightest neutral boson after the photon (Z 1 ) and the electroweak Z boson (Z 2 ). In the M4DCHM 5−5−5 , the lightest composite U and D resonances share approximately the same mass range, from 1.2 TeV to 5.3 TeV, with U 4 typically being the lighter of the two. Similar lower bounds apply to the next three lightest of these species, as well as the lightest of each exotic species, all of which mostly lie below ∼6 TeV but can range upwards of 10 TeV. The remaining fermions can be as light as 2 TeV or 3 TeV, but often lie well above 10 TeV. There are no apparent correlations between the masses across fermion species, aside from the ones that follow from Figure 2. For all valid points, Z 3 is the lightest composite gauge boson, ranging from 1.4 TeV upwards, while Z 4 is approximately degenerate with W 2 and W 3 , which lie above 2.4 TeV. All other gauge bosons tend to be heavier than 10 TeV.
For the M4DCHM 14−14−10 , recall from Figure 2 that two pairs each of U, D, and Q 5/3 particles, and one pair each of Q 4/3 and Q 8/3 particles, come with masses of approximately M ± (m q , m t , m Yt , 0). These are found to typically be among the lightest composite resonances, with masses respectively lying above 1.4 TeV and 1.8 TeV, so there is much approximate degeneracy among the lightest states. This is excepting the U 4 resonance, 9 Not including the gluon sector.  Such resonances would be discovered in collider experiments through the observation of excesses in the invariant mass distributions of the resonance decay products. For the valid points we have found, we give the cross sections of such processes involving the lightest composite resonances that would be measured at the 13 TeV LHC in Figures 12 to 15.
It is seen from Figure 12 that the majority of valid points have predicted U 4 and D 4 cross sections that are far below current collider sensitivity. Exceptions exist in the U 4 → bW + decay channel for the M4DCHM 5−5−5 and the M4DCHM 14−1−10 , and in the D 4 → bZ channel for the M4DCHM 5−5−5 , where even a slight increase in sensitivity would be sufficient to rule out resonances up to ∼1.5 TeV. The exotic decay channels, in Figure 13, are more promising, with all three models making the same remarkably precise predictions for the Q 4/3 1 →bW + and Q 5/3 1 → tW + decays, with cross sections comparable to the upper bounds towards the lighter end of the mass spectrum. Probing these channels further can easily confirm or rule out exotic resonances of certain masses and thereby constrain the models rather straightforwardly, although the cross sections fall off rapidly at higher masses, so much of the space remains out of reach of the LHC. Based on Figure 14, we would not expect to detect any charged composite gauge boson decays in the analysed channels under any of the models. We also considered W 2 → tb decays, but the predicted cross sections for this process were all negligible and not worth presenting. There is better hope for the neutral gauge boson decays in Figure 15, with all models predicting similar cross sections. Promising channels here are Z 3 → tt, and to a lesser extent Z 3 → W + W − . However, many points predict resonance masses too large to be probed in the near future, so failure to detect any resonances at the LHC, for example, would not be sufficient to rule out any of these models. See Refs. [54,55] for more in-depth analyses of the detection prospects for vector bosons in the M4DCHM 5−5−5 .

Higgs signal strengths
Other important theoretical predictions of the M4DCHMs that we analyse are the (gluonfusion produced) Higgs signal strengths. Modifications to the Higgs couplings in a composite Higgs framework have been subject to considerable analysis. To summarise the points that concern us here, in CHMs there are three sources of new physics that can modify the gluon-fusion production of a Higgs boson compared to the SM: non-linearities of the pNGB Higgs (dependent only on v/f ), modified Yukawa couplings of SM particles, and loop contributions from composite resonances [9]. If the contributions from the Higgs nonlinearities were to dominate, as might be expected in our models [56], the signal strengths would be of the form for some channel-dependent (and model-dependent) constants C X [9]. Further discussion on this analytic approximation and the subleading contributions to the signal strengths can be found in Ref. [26]. We give the signal strengths 10 of the valid points in each of our models in Figure 16. The signal strengths µ gg ZZ and µ gg W W are approximately equal due to the custodial symmetry of the M4DCHMs. It is immediately apparent in Figure 16 that the approximation Equation (5.1) does not hold well, especially in the M4DCHM 5−5−5 , but there do seem to be some structures of points that loosely follow the −C X v 2 /f 2 deviation from unity, which is most clearly seen for µ gg γγ in the M4DCHM 14−14−10 and M4DCHM 14−1−10 . This is in contrast to Ref. [27], where clear curves of points were found. It is possible that our points do not follow this distribution cleanly because we consider larger values of f , so the curves are less pronounced and deviations from Equation (5.1) are more apparent.
Among the valid points in the M4DCHM 5−5−5 , the predicted signal strengths in Figure 16 span a wide range of values from 0.8 to 1.3. This agrees with the results of Ref. [57], where it was found that the enhancement in this channel is primarily due to the reduction of the Higgs total width from down-type quark mixing effects. It is therefore unlikely that further probing the signal strengths could rule out or provide evidence for the M4DCHM 5−5−5 . Contrastingly, all the valid points in the M4DCHM 14−14−10 and M4DCHM 14−1−10 give signal strengths that are less than the SM prediction of unity, excepting a small portion of points for µ gg W W and µ gg ZZ . Furthermore, there is a remarkable prediction in these two models: that the signal strength µ gg γγ lies in a narrow range from ∼0.8 to ∼0.9. These 10 The signal strengths for the h → τ τ decay channel are not presented because this channel is subject to large uncertainties, but the predicted values tend to be similar to those for the h → γγ channel.
-38 - clear predictions may serve as tests of the M4DCHM 14−14−10 and M4DCHM 14−1−10 upon more precise measurements of Higgs decays. It is projected that µ gg γγ will be measured to an uncertainty of around 5% by each ATLAS and CMS once they achieve integrated luminosities of 3000 fb −1 in the future high-luminosity run of the LHC [58,59]. At the same time, µ gg ZZ will see a similar improvement in its measurement, while µ gg W W will have only modest improvements. So the diphoton decay channel will serve as the best test of the M4DCHMs, and it is highly likely, based on current measurements, that this channel will disfavour both of these models in the coming years.
It should be stressed that these are the predictions we have found by scanning over the most natural regions of parameter space, under the definition of naturalness encoded by our priors, where low values of f are favoured. It could be the case that more realistic signal strength predictions occur at higher values of f . If the above predictions were to be found false, then, that would only rule out the most natural regions of these models, and not necessarily the models outright. But again, there is a trade-off between accuracy and naturalness that guides the search of CHMs. It would be interesting to investigate the signal strength predictions at larger f in future work.
Although the M4DCHM 14−14−10 and the M4DCHM 14−1−10 give the same predictions for µ gg γγ , they arrive at these predictions for different reasons. In the case of the M4DCHM 14−14−10 , it is the direct collider search constraints, in combination with the other implicit and explicit constraints, that seemingly require µ gg γγ to be significantly less than unity: more realistic signal strengths µ gg γγ 0.9 only occur for otherwise valid points when there are low-mass resonances that are excluded by experiment, as shown in Figure 17. Such a drastic fall-off of µ gg γγ for larger resonance masses is not seen in the M4DCHM 5−5−5 , allowing it to have more realistic signal strengths. For the M4DCHM 14−1−10 , on the other hand, it is the SM mass constraints that favour the lower values of µ gg γγ -particularly the bottom mass, and the top quark and Higgs masses to a lesser degree. Satisfying all three mass constraints simultaneously, or fixing the correct quark Yukawa couplings and Higgs mass, results in all valid points having µ gg γγ 0.9. The fact that the SM mass constraints were included in the scans, while the direct search constraints were not, might explain why the M4DCHM 14−1−10 had more difficulty satisfying the signal strength constraints than the M4DCHM 14−14−10 . But note that if this unfair advantage were to be fixed, the disparity between the evidences of the M4DCHM 14−14−10 and the M4DCHM 14−1−10 would only be exacerbated.

Conclusions
We have performed extensive numerical explorations of the M4DCHM 5−5−5 , M4DCHM 14−14−10 , and M4DCHM 14−1−10 that constitute the first convergent global fits of realistic CHMs. Our fits constrain the parameter spaces into those regions that best reproduce, collectively, the SM masses, the electroweak scale, electroweak precision observables, and various Z boson decay ratios and Higgs signal strengths. Software developed for Refs. [27,29] was employed to calculate these observables, which we modified with updated experimental values. The fits are carried out under a Bayesian statistical framework, aided by the advanced nested sampling program PolyChord, which uses slice sampling to explore the spaces in search of the most natural viable regions.
We determined the Bayesian evidence for each model and took a novel approach in measuring the fine-tuning within a composite Higgs framework, quantifying it with the Kullback-Leibler divergence between the posterior and the prior. With the imposed priors, the M4DCHM 5−5−5 was found to have a Bayesian evidence many orders of magnitude greater than those of the other models, and both the EWSB scale and the top quark mass are incorporated more naturally. Tension between the Higgs signal strength and SM mass constraints was the main obstruction to a higher evidence for the M4DCHM 14−1−10 , while a considerable fine-tuning resulted in the M4DCHM 14−14−10 having the lowest evidence. We found that these results are strongly prior dependent and we stress that a future conclusive model comparison requires a detailed analysis of the impact of the prior distributions and bounds. Interestingly, it was found that the posterior occupies surfaces of only ∼7 − 12 effective dimensions in the (18+ dimensional) parameter spaces of the models.
Each model was found to be capable of satisfying all constraints individually to within 3σ, and their collider phenomenology was analysed in these viable regions. Our results suggest that the fermionic resonances in these models are excluded below ∼1.1 TeV. We expect the models to have different signals for the various decay modes of the up-type and down-type resonances at the √ s = 13 TeV LHC, but remarkably all models give the same precise predictions for the exotic fermion decays Q 5/3 → tW + and Q 4/3 →bW + . These exotic decays constitute the most promising channels for probing these models in future collider searches. The lightest composite vector bosons range in mass from ∼1 TeV to over 10 TeV, and neutral bosons on the lighter end of the spectrum can be probed effectively in the Z 3 → tt decay channel. Finally, the Higgs signal strengths predicted by -40 -the M4DCHMs were analysed, and the models M4DCHM 14−14−10 and M4DCHM 14−1−10 were found to have a clear signature: they predict the gg → H → γγ cross section to be between ∼80% and ∼90% of the value predicted by the SM. This is already in tension with experiment at around the 2σ level, and has the potential of being strongly ruled out in the future high-luminosity run of the LHC.

Acknowledgements
MJW, AGW and WS are funded by the ARC Discovery Project DP180102209, the ARC Centre of Excellence for Dark Matter Particle Physics CE200100008 and the Centre for the Subatomic Structure of Matter (CSSM).

A SO(5) generators and representations
Here we give the conventions for SO(5) generators and embeddings we use in this work. This section makes extensive use of local isomorphism SO(4) ∼ = SU (2) L × SU (2) R , since this allows representations of SO(5) ∼ = SO(4) × SO(5)/SO(4) to be decomposed into those of SU (2) L × SU (2) R , making the embedded fields' quantum numbers easily identifiable.

Generators
To exploit the local isomorphism, the generators are conveniently split into (T a L,R ) 3 a=1 , which generate the subgroup locally isomorphic to SU (2) L,R , and the broken generators (T b ) 4 b=1 associated with the SO(5)/SO(4) coset. In the fundamental representation they are given by Note this is consistent with our choice of the SO(4)-invariant vector Φ 0 = (0, 0, 0, 0, 1) that implicitly embeds SO(4) into the top left 4×4 block of SO(5) matrices.

Field Embeddings
The representations of SO(5) that we are interested in decompose under SO(4) and subsequently under SU (2) L × SU (2) R as Multiplets in the representation r will be denoted Ψ r . Representations of SO(4) are embedded inside those of SO(5) as [60] Here we are treating Ψ 10,14 as matrices, which are acted upon by the generators through the commutator: The SO(4) multiplets are populated with fields Ψ n L ,n R having SU (2) L × SU (2) R quantum numbers 11 (n L , n R ) according to where Ψ 6 is antisymmetric and Ψ 9 symmetric. The superscript ± denotes ±1/2 in Ψ 4 and ±1 in the other multiplets, and we have defined Two different fields have (n L , n R ) = (0, 0) in Ψ 6 .

B Correlators
Here we provide expressions for the correlators Π t,b L,R and M t,b for each model of the M4DCHM, originally given in Ref. [26]. See Equation ( In this model, The form factors for the bottom sector are obtained by interchanging all t andt subscripts with b andb.

C Mass matrices
Here we give the explicit mass matrices for each model of the M4DCHM. They are given in the Site 0 holographic gauge, in bases consisting of the fields with definite SU (2) L ×SU (2) R quantum numbers presented in Appendix A. Only particles of the same electric charge can mix, so the mass matrices are separated according to the particles' charges.

C.1 Boson sector
The mixing terms between the gauge bosons are independent of the fermion representations, so they are the same across all models considered in this work. Their mixings are represented by symmetric matrices, whose singular values are the squared masses of the resonances. Expansion of the Lagrangian in Equation (2.3) reveals the SU (3) resonances only mix among themselves, through the matrix Note one singular value of this matrix is 0, as expected for the massless SM gluons. Of the remaining bosons, four are charged and seven are uncharged. The charged bosons {W 0 ± µ , ρ ± Lµ , ρ ± Rµ , a ± µ } are linear combinations of the gauge fields that come paired with the generators (using our choice of generators) given by the form Their mass matrix is given by Equation (C.3), and the uncharged bosons have the mass matrix Equation (C.4).

C.2 Fermion sector
The Lagrangians Equations (2.19), (2.22) and (2.24) lead to the quark partner mass matrices for each model below. Subscripts denote the matrices for the up-type (U) and down-type (D) fields, as well as the fields of exotic charges 4/3, 5/3, and 8/3. The subscripts/tildes of the fields match those of the multiplets in which they are contained. Where necessary to avoid ambiguity, subscripts denote the representations of SO (5) or SO(4) in which the fields reside. Fields may share a label with those in other matrices due to the different definitions of ± in superscripts for different SO(4) representations (see discussion after Equation (A.8)), but the electric charges distinguish them and remove ambiguity about the representation each comes from. To save space in writing the matrices, we define the quantities and A caveat: in this model a different convention was used for the field embeddings than was presented in Appendix A. The fields in Ψ 4 in Appendix A have been redefined here as Ψ −,− → −iΨ −,− , Ψ +,+ → iΨ +,+ , and Ψ +,− → −iΨ +,− . (C.10) -48 - (C.14) -50 -       Our results for this model come from four different PolyChord scans. The first two scans, which each used 2000 live points, did not display very good agreement with each other, but more agreement was seen after increasing the sampling density to 4000 live points. The posteriors found by these 4000 point scans are shown in Figures 18 to 20. There is not exact agreement in these posteriors -for example, the second run found a minor mode that the first did not, apparent in Figures 19 and 20 -but there are no major differences in the general regions that each scan found, and specific features such as the correlation between m b and m Y b are present in both scans. Parameters for which there is poor agreement, such as g X and g G in Figure 18, tend to be those that are are minimally affected by the constraints we impose and as such are not our primary interest. Given how difficult sampling the space has proven to be, we regard this as an acceptable level of agreement. It is also encouraging that the Bayesian evidences found in the 4000 point scans were in agreement: which we take as further indication that the parameter space of this model has been reliably explored.
-57 - As such, these results are quite reliable.