{'fraction_non_alphanumeric': 0.03383045717493715, 'fraction_numerical': 0.018867924528301886, 'mean_word_length': 4.617715617715618, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Bayesian rotation inversion of KIC 11145123

A scheme of Bayesian rotation inversion, which allows us to compute the probability of a model of a stellar rotational profile, is developed. The validation of the scheme with simple rotational profiles and the corresponding sets of artificially generated rotational shifts has been successfully carried out, and we can correctly distinguish the (right) rotational model, prepared beforehand for generating the artificial rotational shifts, with the other (wrong) rotational model. The Bayesian scheme is applied to a gamma Dor-delta Sct type hybrid star, KIC 11145123, leading to a result that the convective core of the star might be rotating much faster (~ 10 times faster) than the other regions of the star. The result is consistent with that previously suggested by Hatta et al. (2019) based on a 3-zone modeling, further strengthening their argument from a Bayesian point of view.


INTRODUCTION
Stellar internal rotation plays an essential role in a variety of stellar internal dynamics (e.g. Maeder 2009) such as generation of magnetic fields via the dynamo mechanism and transportation of chemical elements caused by rotationally induced mixing. It is therefore of great importance for us to study stellar internal rotation theoretically and observationally, the latter of which, in particular, has recently become feasible thanks to the establishment of asteroseismology (e.g. Aerts et al. 2010) brought about by high-precision photometric observations by spacecrafts such as Kepler (Koch et al. 2010) and TESS (Ricker et al. 2014), leading to numerous asteroseismic inferences on internal rotation of various types of stars including solar-like stars (e.g. Benomar et al. 2015Benomar et al. , 2018Schunker et al. 2016a,b), early-type main-sequence stars (e.g. Kurtz et al. 2014;Saio et al. 2015;Schmid and Aerts 2016;Papics et al. 2017;Christophe et al. 2018;Ouazzani et al. 2019;Li et al. 2020), and evolved stars (e.g. Beck et al. 2012;Mosser et al. 2012;Di Mauro et al. 2018;Deheuvels et al. 2020).
The current understanding of internal rotation of main-sequence stars is summarized by Aerts et al. (2019) that almost all the main-sequence stars investigated so far are exhibiting nearly rigid rotation throughout them, which has not been expected based on previous hydrodynamical numerical simulations of angular momentum transfer inside stars (e.g. Tayar and Pinsonneault 2013;Eggenberger et al. 2017).
To fill the gap between the observation and the theory, several mechanisms of angular momentum transfer by, for instance, internal gravity waves or magnetic fields, have been proposed (e.g. Cantiello et al. 2014;Rogers 2015;Fuller et al. 2019). Thus, asteroseismic studies have definitely been contributing in propelling understanding of the stellar rotation.
Interestingly, there are also a few asteroseismic researches suggesting the existence of rotational velocity gradient inside stars (e.g. Benomar et al. 2018); stars are rotating rigidly throughout most of the interiors, but not completely. One of such stars for which rotational velocity shear inside has been suggested is KIC 11145123 (Hatta et al. 2019), which is a γ Dor-δ Sct type hybrid star (Bradley et al. 2015) and has been actively studied based on its well-resolved frequency splitting for pressure (p), gravity (g), and mixed modes, revealing the evolutionary stage (Kurtz et al. 2014;Hatta et al. 2021), the asphericity (Gizon et al. 2016), and the internal rotation (Kurtz et al. 2014;Hatta et al. 2019). Hatta et al. (2019)'s primary focus was on investigating the latitudinally differential rotation of the star, but, by carefully checking behaviors of estimates obtained via rotation inversion, they have found a hint that the convective core of the star might be rotating 5-6 times faster than the other regions of the star. Although the suggestion of the fast-core rotation might appear to be incompatible with the current understanding that mainsequence stars are rotating almost rigidly, it is actually not the case. Previous asteroseismic studies, especially those focusing on early-type main-sequence stars with the convective core and the broad radiative region above, have utilized high-order g modes, which are not established in the convective core, to infer rotation rates in the deep region; what they have estimated are rotation rates in the deep radiative regions. In contrast, Hatta et al. (2019) have used mixed modes, which have finite sensitivity, though small, inside the convective core, enabling them to obtain a hint of the rotation rate of the convective core located beneath the radiative region. They have actually confirmed that the star is rotating almost rigidly throughout the radiative region and that the convective core is the exception. Saio et al. (2021) provide another asteroseismic study which have inferred convective-core rotation of fastrotating γ Dor stars by fitting characteristic dips in observed g-mode period spacings, caused by the coupling between inertia modes and high-order g modes. Though they have found no hint of rotational velocity shears (between the convective core and the radiative region above) among their targets, comparison between the two studies should be helpful for putting further constraints on theoretical calculations of angular momentum transfer inside stars.
In this paper, we would like to take a further step in terms of investigating convective-core rotation of stars based on a newly developed scheme of Bayesian rotation inversion which enables us to compute probabilities of models of rotational profiles based on the so-called global likelihood (e.g. Gregory 2005). The first goal is to present the scheme of Bayesian rotation inversion, and the second goal is to apply the scheme for KIC 11145123 and to examine if the fast-convective-core rotation is obtained or not.
The structure of the paper is as follows. In Section 2, the mathematical formulations of the Bayesian rotation inversion is presented after a brief introduction to asteroseismic rotation inversion, which is based on the perturbative approach. Note that the perturbative approach is justified for inferring the internal rotation of KIC 11145123 because of the fact that the star is a slowrotator with the rotation period of about 100 d, which is much longer than the dynamical timescale as well as the oscillation periods of the star of a few hours. In Section 3, validation of the developed scheme is carried out with simple artificial datasets. We apply the Bayesian scheme to KIC 11145123 in Section 4 where the basic properties of the star, the comparison of modeled rotational profiles, and results finally obtained have been featured. We lastly give a summary in Section 5.
Finally, we have a note on reference stellar models used in this paper. As reference stellar models for rotation inversion, we firstly chose two models; one is a model constructed by Kurtz et al. (2014) assuming single-star evolution, which was used in Hatta et al. (2019), and the other is a non-standard model of the star constructed by Hatta et al. (2021) taking effects of some interactions with other stars into account. We, however, did not see any significant differences in results of rotation inversion no matter which model was taken as a reference model, and thus, we will present the results obtained based on the newer model computed by Hatta et al. (2021) in this study. Some of global stellar parameters of the reference model are as follows: M = 1.36M ⊙ , Y init = 0.26, Z init = 0.002, Age = 2.160 × 10 9 years; the star is represented by a relatively low-mass stellar model around the terminal age main-sequence stage, exhausting most of the hydrogen at the central hydrogen-burning region. For more information on the reference model, see Hatta et al. (2021).

Rotation inversion
It is relatively simple to mathematically describe the frequency splitting caused by rotation (the rotational splitting) when the internal rotation is slow compared with the dynamical timescale and the oscillation periods of the system. In this case, we can treat rota-tion as a small perturbation to the system, and then, we can relate the internal rotation to the rotational shifts in frequencies based on the first-order perturbative approach (see more details in e.g. Unno et al. 1989;Aerts et al. 2010). The explicit form for the rotational shift (d i = δω i /m, where i represents a particular set of mode indices, namely, the radial order n, the spherical degree l, and the azimuthal order m) thus derived is as follows: in which the internal rotation Ω(x, µ) is expressed as a function of a position inside the star, represented by the fractional radius (x = r/R) and the cosine of the colatitude θ (µ = cos θ). The observational uncertainty for the rotational shift is given by e i . The rotational splitting kernel K i (x, µ) can be obtained by calculating the linear adiabatic oscillation of a certain reference model, and the explicit form can be found in e.g. Aerts et al. (2010).
Then, what we have to do to estimate the internal rotation Ω(x, µ) is to invert the set of the equations (1) (rotation inversion), where the number of the equations is identical to that of observed rotational shifts. Techniques such as the Regularized Least-Squares (RLS) method (e.g. Tikhonov and Arsenin 1977) and the Optimally Localized Averaging (OLA) method (Backus and Gilbert 1967), both of which are well established, have been frequently used in helioseismology, which is also the case in asteroseismology.
Note that, in asteroseismology, we do not have a large number of rotational shifts (around a few dozens at most) compared with the case in helioseismology (of the order of 10 5 ), and thus, it is sometimes difficult to draw definitive conclusions based on just one method (even when it is one of the standard methods) and comparisons of results obtained via different methods can help us to better understand the inversion results, which is another reason why we attempt to develop a new scheme of Bayesian rotation inversion in this study. In other words, there is no all-round method which enables us to solve any inverse problems completely, i.e., each inversion technique provides us with the corresponding estimate based on a particular criterion adopted for the technique. Therefore, we should be cautious not to jump to seemingly satisfactory conclusions, which is especially the case in asteroseismology where the relative scarcity of observed rotational shifts easily leads to the ill-posedness of the inverse problems.

Basic points in Bayesian statistics
In this section, we would like to give a few basic points in Bayesian statistics. Bayesian statistics allows us to, for example, investigate global properties of probabilities of parameters or compute probabilities of models, based on the latter of which we can conduct model comparison among possible models. In particular, such capability of the model comparison strongly motivates us to construct an inversion scheme based on Bayesian statistics, and the application can be found in Section 4, with which the possibility of the fast-convective-core rotation of KIC 11145123 is tested. For more thorough introductions and discussions on Bayesian statistics in astronomical contexts, readers should refer to e.g. Gregory (2005). Applications in global and local helioseismology can be also found in Kashyap et al. (2021) and Jackiewicz (2020), respectively.
Then, let us introduce one of the most fundamental equations in Bayesian statistics, the Bayes' theorem, which has the following form: which can be derived based on the definition of the conditional probability. Datasets (obtained by observation) and parameters (to be estimated) are represented by d and θ following the notation in Benomar et al. (2009). The probabilities p(θ|d), p(d|θ), p(θ), and p(d) should read the posterior probability of the parameters given the dataset, the probability of the dataset given the parameters (or, the likelihood of the parameters), the prior probability of the parameters, and the probability of the dataset marginalized by all the parameters (the so-called global likelihood), respectively. What the Bayes' theorem (2) indicates is actually not complex to interpret; though we have to begin with uninformative states (represented by the prior probability p(θ)), once we conduct observations, we can estimate the probability of obtaining the resultant datasets assuming a set of parameters (represented by the likelihood of the parameters p(d|θ)), and finally, combining the prior probability and the likelihood enables us to update our understanding of the parameters (represented by the posterior probability p(θ|d)).

Formulation
One of the goals in the Bayesian rotation inversion is to compute the posterior probability of parameters describing a rotational profile given the observed rotational shifts p(Ω(x, µ)|d). To this end, we need to determine functional forms of the prior probability p(Ω(x, µ)) and the likelihood p(d|Ω(x, µ)) beforehand as shown in Subsection 2.3.1. We also describe how to carry out model comparison in the Bayesian framework in Subsection 2.3.2. In this subsection, we basically follow the formulations in Benomar et al. (2009).

Priors and likelihoods
For the prior probability, we adopt a uniform distribution whose range is reasonably broad, as we, a priori, do not have strong constraints on the parameters to be estimated. Although it is generally recognized that the choice of functional forms of prior probabilities could affect Bayesian inferences (e.g. Benomar et al. 2009), this seems not to be the case in our study as will be seen in Section 4, which is another reason for adopting simple uninformative prior probabilities for the parameters. The specific range for the prior probability is later shown in Sections 3.2 and 4.2.
For the likelihood, we assume that an observed rotational shift d i is a realization from a Gaussian distribution whose mean is given by the first term on the righthand side in expression (1) with its standard deviation identical to that of the corresponding observational uncertainty e i . Based on the assumption, the explicit form for the likelihood of the parameters given the rotational shift is as below: By further assuming that the observed rotational shift is statistically independent of each other, the explicit form for the likelihood of the parameters given the set of the observed rotational shifts is the product of expressions (3) as follows: where the number of the observed rotational shifts is denoted by N . Based on the determined prior probability and the likelihood of the parameters given the set of the observed rotational shifts, we can calculate the posterior probability of the parameters following expression (2), and we can subsequently obtain estimates for the parameters by, for instance, choosing a set of the parameters for which the posterior probability is maximum (called Maximum A Posteriori estimation).
Note that we have not determined an explicit functional form for the rotational profile Ω(x, µ) yet, which is a necessary step for us to compute the likelihood (we have to compute the integration inside expression 3).
Several specific rotational profiles (and the corresponding results of the Bayesian rotation inversion) can be found in Sections 3 and 4. It should also be noted that calculation of the posterior probability requires us to carry out numerical integrations via, for instance, the Markov Chain Monte Carlo (MCMC) method if the number of the parameters used to describe a rotational profile is so large that it is computationally impossible to directly evaluate the posterior probability.

Model comparison based on global likelihoods
In this subsection, we would like to mention a model comparison based on Bayesian statistics. The important quantity is the global likelihood p(d), which is a normalization constant in expression (2). We can confirm the importance of the global likelihood by reconsidering the Bayes' theorem (2), which can be rewritten as where a model M j representing a certain set of parameters is explicitly expressed, and the global likelihood can read the likelihood of the model M j given the dataset.
Then, let us consider the posterior probability of the model M j as below: in which the ratio of the posterior probabilities of the models is expressed by the ratio of the global likelihoods of the models (note that it is assumed that Jeffreys 1998). This threshold can be explained with a simple example where there are just two models M j and M k . In that case, the posterior probability of the model M j can be expressed as Mj ,M k ] −1 . Thus, the threshold 3 corresponds to the probability 0.75. In this way, the global likelihoods p(d|M ) are such essential quantities that we can select the most probable model given the dataset, which is practically demonstrated in Sections 3 and 4.

SIMPLE TESTS
The Bayesian scheme explained in Section 2 is tested in this section. We firstly introduce two models of rotational profile, from which two specific rotational profiles, one with a rotational velocity shear and the other without it, are artificially generated. The corresponding sets of rotational shifts are given as well (Section 3.1). We carry out Bayesian rotation inversion with the artificial datasets to compute the posterior probabilities and global likelihoods, based on which whether we can choose the correct model of rotational profile or not is checked (Section 3.2).

Artificial rotational shifts
A rotational shift of a certain mode d i can be computed following expression (1) once we specify a "true" rotational profile Ω(x, µ) and calculate the rotational splitting kernel K i (x, µ) of the model. We present two models of rotational profile Ω(x) as below (note that we do not consider latitudinal dependence of the internal rotation for simplicity in this section). The first one is described with two parameters, namely, a rotation rate of the core Ω 1 and that of the envelope Ω 2 (the red lines in Figure 1). The other one has a linear profile in terms of the fractional radius, which is parameterized with a rotation rate of the center of the star Ω c and the surface rotation rate Ω s (the blue line in Figure 1). Let us call the former model M sh and the latter one M lin . Specific parameters for two particular rotational profiles are here given as below: (Ω 1,true , Ω 2,true ) = (1.50, 0.980), and (Ω c,true , Ω s,true ) = (0.950, 0.980), in units of 2π × 0.01 d −1 . Let us call the rotational profiles thus specified Ω sh and Ω lin , respectively. Note that, in this study, we concentrate on mainsequence stars with a convective core as is the case for KIC 11145123. Accordingly, we assume that the boundary in Ω sh is fixed to be the convective boundary; Ω 1 and Ω 2 represent the rotation rate of the convective core and that of the radiative envelope.
For computing splitting kernels K i (x, µ), the linear adiabatic oscillation of the reference model, which is the non-standard model of KIC 11145123 constructed by Hatta et al. (2021) (see Section 1), has been calculated via a linear adiabatic oscillation code GYRE (Townsend and Teitler 2013). Based on the eigenfunctions and eigenfrequencies thus obtained, we calculate the splitting kernels K i (x, µ) (the explicit form can be found in, e.g., Aerts et al. 2010). We have computed splitting kernels for 33 eigenmodes, namely, 20 profile Ω(x), namely, M sh (red) and M lin (blue), which are described with the sets of two parameters (Ω 1 , Ω 2 ) and (Ωc, Ωs), respectively.
high-order g modes with (l, m) = (1, 1), 3 low-order p modes with (l, m) = (1, 1), 5 low-order mixed modes with (l, m) = (2, 1), and 5 low-order mixed modes with (l, m) = (2, 2). The eigenmodes are thus chosen since these types of modes (high-order g modes and low-order p/mixed modes) are frequently observed for γ Dor-δ Sct type hybrid stars such as KIC 11145123 though the relatively larger number of modes compared with the actual observation are prepared to render rotation inversion as robust as possible. It should be noted that the results of the simple test here are qualitatively the same even if we use a smaller set of splitting kernels, which is identical to that in the case of KIC 11145123 (see Section 4.1).
Then, based on the rotational profiles, namely, Ω sh and Ω lin , and the splitting kernels K i (x, µ), we have computed the corresponding sets of rotational shifts following expression (1), in which the observational uncertainties e i are assumed to be realizations from a Gaussian distribution whose mean and standard deviation are 0 and 10 −3 (in units of 2π × 0.01 d −1 ). The standard deviation of the Gaussian distribution is chosen so that it is around a typical observational uncertainty of the rotational shifts for KIC 11145123. It is also assumed that the observational uncertainty e i is statistically independent of each other. We thus have two sets of 33 artificially generated rotational shifts with which Bayesian rotation inversion is to be carried out in the following sections. We denote the artificial dataset for the rotational profile with a velocity shear as δω sh and that for the linear rotational profile as δω lin .
3.2. Bayesian rotation inversion with the artificial rotational shifts  In this section, the Bayesian scheme has been utilized with the set of rotational shifts (δω sh or δω lin ), in order to test whether we can correctly choose the right model of rotational profile (M sh or M lin , respectively) via the Bayesian scheme or not.
As mentioned in Subsection 2.3.1, the first thing we have to do is to specify parameters to describe Ω(x, µ) for which the posterior probabilities are to be computed. For simplicity, we have used the models of rotational profile M sh described by (Ω 1 , Ω 2 ) and M lin described by (Ω c , Ω s ) (Figure 1), the former of which contains Ω sh and the latter of which contains Ω lin .
Secondly, prior probabilities for the prepared parameters have to be specified. We assume that the prior probabilities are uniform as below: where q ∼ U [a, b] means that q is a random variable uniformly distributed in a range from a to b. The rotation rates are in units of 2π × 0.01 d −1 . The joint prior probability is computed by taking products of prior probabilities assuming that the parameters are statistically independent of each other. Based on the priors and likelihoods, which can be computed with expressions (3) and (4), we calculate the numerator of the right hand side in expression (2)   Note. Based on the definition (7), it is apparent that It is generally considered that the model M i is substantially favored compared with the model M j when log O M i ,M j > 0.5. malization constant p(d) (or, the global likelihood), and finally compute the posterior probability of parameters given the dataset. Since the numbers of the parameters are just two for both of the models and it is not computationally expensive to numerically carry out such 2dimensional computations, we have directly calculated the posterior probabilities. Figure 2 shows two posterior probability density functions given δω sh . One is computed based on M sh and the other is computed based on M lin (P sh and P lin , respectively, in Figure 2). Since both of the posterior probabilities are unimodal, we can estimate the parameters by, for instance, taking a parameter set that maximizes the corresponding posterior probability (the Maximum A Posteriori estimate). Figure 2 that the Maximum A Posteriori estimate for P sh is almost identical to the prepared parameters (8) (red dotted lines in Figure 2) while those for P lin are biased. The important point, however, is that we cannot determine which model describes the dataset δω sh more appropriately by just comparing the posterior probabilities without knowing the prepared parameters (8) beforehand. As it is described in Subsection 2.3.2, such model comparison should be achieved by comparing the global likelihood, which is the normalization constant in expression (2).

It is apparent in
In the case of this simple test with δω sh , we have obtained the following odds ratio (defined by equation (7)): log O M sh ,M lin = 23.10 (Table 1), which is much larger than the conventionally adopted criterion log 3 ∼ 0.5 (e.g. Jeffreys 1998). It is therefore correctly inferred that, in the light of the Bayesian scheme, the model of rotational profile with a velocity shear M sh is more favorable to describe the artificially generated dataset δω sh .
The same is true when δω lin is used for the test (Figure 3). The obtained odds ratio is log O M lin ,M sh = 135.6( Table 1). It is clear that the model M lin is preferred, which is the right one we have used to generate the dataset δω lin . It should also be instructive to mention that the estimates are biased unless we have chosen the correct model (see P sh in Figure 3), highlighting the importance of the model comparison in the Bayesian context achieved by computing global likelihoods.

APPLYING THE METHOD TO KIC 11145123
The Bayesian scheme demonstrated in the previous sections is applied to one of the Kepler targets, KIC 11145123, to infer its internal rotation profile, especially focusing on the convective-core rotation. After we present the rotational shifts and corresponding splitting kernels for KIC 11145123 (Section 4.1), basic setups for the Bayesian rotation inversion are given (Section 4.2), based on which the posterior probabilities and the global likelihoods are computed, and the model comparison is conducted as well (Section 4.3). Then, the validation of the results obtained is carried out in Section 4.4. We finally provide a brief discussion on the results in Section 4.5.

Data and splitting kernels
We have used a set of rotational shifts and observational uncertainties which have been measured and determined by Kurtz et al. (2014) (see Tables 1 and 2 in the paper). The set is composed of twenty-three eigenmodes, namely, fifteen high-order g modes with (l, m) = (1, 1), two low-order p modes with (l, m) = (1, 1), three low-order mixed modes with (l, m) = (2, 1), and three low-order mixed modes with (l, m) = (2, 2). The mode identification for the mode set has been conducted based on the non-standard model of the star constructed by Hatta et al. (2021), according to which the star is a low-mass star at the terminal age main-sequence stage. Some of basic global parameters of the non-standard model can also be found in the last paragraph of Section 1. The corresponding splitting kernels are calculated using the non-standard model as a reference model. We have three types of the splitting kernels. The first type corresponds to high-order g modes and it has sensitivity in the deep radiative region just above the convective core (Figure 4). The second type and the third correspond to low-order p modes and low-order mixed modes, respectively (Figures 5 and 6). Both of them have sensitivity in the outer envelope, but the low-order mixed modes with (l, m) = (2, 1) have sensitivity in the high-latitude region as well, which is not the case for other modes with (l, m) = (1, 1) or (l, m) = (2, 2). It should be emphasized that only the mixed-mode splitting kernels have finite sensitivity inside the convective core, which renders the detection of the convective-core rotation possible as it has been shown in Hatta et al. (2019). A close look into the splitting kernels around the convective boundary is later presented and discussed in Section 4.5. profile as below:

Models of rotational profile and priors
where Ω 0 (x) can have a linear profile described with two parameters, namely, the rotation rate at the center Ω c and that at the surface Ω s (see the blue line in Figure  1) or can have a velocity shear (with four parameters, namely, the rotation rate below a velocity shear boundary x sh which is assumed to be uniform at Ω core , the rotation rate Ω rad at another side of the shear boundary fixed to be the convective core boundary (x czb ∼ 0.045), and that at the surface Ω s ) (see Figure 7). Another function Ω 1 (x), which is related to the latitudinal dependence of the internal rotation, can be zero everywhere (with no parameters) or can have a linear profile (with two parameters, in almost the same way as a linear profile of Ω 0 (x) but with additional indices as Ω 1c and Ω 1s ). We have taken such latitudinal dependence of the internal rotation into account in order to evaluate the effect on the global likelihood though it is not a primary subject to be investigated in this study. We note that the final inference on the convective-core rotation is not qualitatively changed due to the inclusion of the latitudinal dependence of the rotation in the analysis as we see later in Section 4.3.
With the definitions for Ω 0 (x) and Ω 1 (x), there are four ways of parameterization of the rotational profile in total. Let us denote the models as follows: M 1d2p for linear Ω 0 and zero Ω 1 , M 1d4p for shear Ω 0 and zero Ω 1 , M 2d4p for linear Ω 0 and linear Ω 1 , and M 2d6p for shear Ω 0 and linear Ω 1 , where the subscripts d and p stand for the dimension of a modeled rotational profile and the number of parameters in the model, respectively (Table  2). Then, what we have to specify is the corresponding priors for the parameters. As described in Section 3.2, each prior probability for a certain parameter is assumed to be uniform as below: and x sh ∼ U [0.010, 0.055].
The parameters representing rotation rates are in units of 2π × 0.01 d −1 , and q ∼ U [a, b] means that q is a random variable uniformly distributed in a range from a to b, as described in Section 3.2. The joint prior probability is computed by taking products of the prior probabilities based on the assumption that the parameters are statistically independent of each other.

Results
For each way of parameterization, the posterior probability of the parameters is computed based on the likelihood of the parameters and the joint prior probability. The likelihood is calculated with expressions (3) and (4). We carry out the so-called Metropolis method (Metropolis et al. 1953), which is one of the standard algorithms to carry out MCMC, to evaluate the posterior probability. The convergence of samples generated via the Metropolis method has been checked by visual inspection (see Figure 8) with several different sets of initial values for the samples. A typical number of iterations required is of the order of 10 5 , and small fractions of samples are discarded from the final samples Ω core Ω rad Ω s x sh x czb Ω x Figure 7. Schematic picture for a way of parameterizing the function Ω 0 (x) in expression (11), which is used for the models M 1d4p and M 2d6p . The profile contains a velocity shear, and it is described by four parameters, namely, the rotation rate inside the convective core Ωcore, that at the outer boundary of the velocity shear Ω rad , that at the surface Ωs, and the position of one of the boundaries of the velocity shear x sh . Note that another boundary of the velocity shear is fixed to be the convective boundary (represented as x czb in this figure); x sh could be either the inside or the outside of the convective core. Except for the inner core which is assumed to be rotating rigidly with Ωcore, rotation rates of the other regions are linearly expressed.
as they are considered as samples in the burn-in periods (a period during which obtained samples are not thought to be realizations from the posterior probability distribution we would like to sample). More information on principles of MCMC, how to manage the outcomes of MCMC, and so on, can be found, for example, in Gregory (2005). We consider that the posterior probabilities have been successfully sampled via MCMC, based on the convergence of MCMC samples (see Figure 8). Then, we compute the global likelihood for each way of parameterization of the rotational profile to carry out the Bayesian model comparison. Note that the global likelihood p(d), which is the normalization constant in relation (2), cannot be determined with the MCMC samples alone. We therefore need some tools to compute the global likelihoods. This can be accomplished by following the procedures proposed by Chib and Jeliazkov (2001) in which an exact value of the posterior probability of a particular set of the parameters is directly evaluated, and then, the global likelihood is calculated based on relation (2). We have confirmed that the method correctly works for simple cases where we can analytically compute the posterior probability and the global likelihood.
The decimal logarithm of the odds ratios thus ob- likelihoods has been checked with several different settings in MCMC, and the global likelihoods have been rounded to no more than four significant figures (Table  2).
We can also confirm that based on the posterior probability for the models with the velocity shear, the fastcore rotation is again inferred. Figure 9 shows the 2dimensional contour map of the (marginalized) posterior probability density function as a function of (Ω conv , x sh ), which is computed for the model M 2d6p . Note that Ω conv is used instead of Ω core to describe the latitudinally averaged rotation rate of the convective core. When we adopt the Maximum A Posteriori estimates for the parameters, they are (Ω conv , x sh ) ∼ (13.4, 0.0493), and thus, the fast-convective-core rotation has been inferred as in the case of Hatta et al. (2019) (see Figure 8 for Ω s and Ω rad , which are of the order of 1 in units of 2π × 0.01 d −1 and are much slower than Ω conv ).

Validity of the results
Since the models used for Bayesian model comparison in the last section are much more complex than those in Section 3, we would like to check whether we can conduct Bayesian model comparison even with such complex models, which eventually enables us to validate the results obtained in the last section. To this end, we present the same test as that described in Section 3 except for the following two points. Firstly, the Maximum A Posteriori estimates for the parameters, which are determined with the posterior probability given the model M 2d6p (computed in the last section), are used for constructing an artificial rotational profile and the corresponding rotational shifts. Let us call the artificial rotational shifts δω 2d6p . Secondly, to compute the posterior probabilities given the artificial rotational shifts, we have used the four models, namely, M 1d2p , M 1d4p , M 2d4p , and M 2d6p (see the definitions of the models in Section 4.2). What we would like to check is whether we can correctly choose the right model (M 2d6p in this case) among the models based on the global likelihoods or not.
The global likelihoods have been computed via the method of Chib and Jeliazkov (2001 (Table 3). It is therefore evident that the correct model M 2d6p is preferred to the other models from the Bayesian perspective. This can be also confirmed when we see Figure 10, where the prepared parameters (represented by red dotted lines in Note. Based on the definition (7), it is apparent that

Discussion
Based on the results obtained in Section 4.3, it has been inferred, via the Bayesian scheme, that the convective core of KIC 11145123 is rotating approximately ∼ 10 times faster than the other regions of the star. We have also confirmed the validity of the results in Section 4.4. In this section, we have a brief discussion about the inferred position of the velocity shear. Figure 11 shows the Brunt-Väisälä frequency of the referencemodel (black curve) and the inferred position of the velocity shear (blue dashed line). What we would like to point out is that the inferred position of the velocity shear x sh ∼ 0.05 is slightly above the convective boundary x czb ∼ 0.045, at which the square of the Brunt-Väisälä frequency is zero, and that the shear is located in the overshoot zone (represented as the blue shaded area in Figure 11). In the overshoot zone, where overhooting has been modeled as a diffusive process following Herwig (2000) in the case of the reference model, the square of the Brunt-Väisälä frequency is positive so that the overshoot zone is included in the g-mode cavity. This is fairly relevant to the inversion analysis since not only the mixed-modes but also the high-order g modes should also have sensitivity for the fast-core rotation (compare the black curve and the orange curve in Figure 12).
The relevance described in the last paragraph can be actually confirmed by a test, in which the Bayesian rotation inversion is carried out without the mixed-mode rotational shifts. We then have the following odds ratio in the decimal logarithmic scale: log O M 1d4p ,M 1d2p = 20, ( Table 2). Note that the latitudinal dependence is not considered here for simplicity. Although the difference in the global likelihoods are much smaller than those Splitting Kernel X high-order g mode mixed mode Figure 12. Slices of the mixed-mode (black) and high-order g-mode (orange) splitting kernels at the equator (µ = 0). The horizontal axis shows the fractional radius x. The inferred position of the velocity shear is represented by the blue dashed line. The convective core and the overshoot zone are represented by the grey shaded area and the blue shaded area, respectively.
computed with the mixed modes included in the analysis, we can still claim that the model with a velocity shear M 1d4p is more favorable in terms of the global likelihood; the inference of the velocity shear has just become more marginal. The same analysis with the artificial rotational shifts δω sh in Section 3.2 has shown that it becomes rather marginal to distinguish the two models of rotational profile (one with a velocity shear (M sh ) and the other without it (M lin )) if we exclude the mixed modes. The odds ratio expressed in the decimal logarithm is log O M sh ,M lin = 1.25. This is mainly because we artificially locate the position of the shear identical to the convective core boundary and there are no modes which are sensitive to the artificially generated fast-convectivecore rotation other than the mixed modes. However, in the case of KIC 11145123, high-order g modes can be sensitive to the fast-core rotation which slightly invade the radiative region, possibly leading to the marginal preference to the model with a velocity shear M 1d4p despite the absence of the mixed-mode rotational shifts in the analysis.

SUMMARY
With the increasing number of stars for which asteroseismic rotation inversion can be carried out, developing and improving inversion techniques is of great importance for us to render our inferences on the internal rotation of stars as robust and reliable as possible. In this study, we first present a scheme of Bayesian rotation inversion which enables us to compute the probability of a model of rotational profile and thus select the most reliable model among multiple models prepared by us beforehand (model comparison via the global likelihood). We then conduct a simple test for the scheme using two models of rotational profile, based on which two specific sets of rotational profiles and the corresponding rotational shifts are artificially generated. It has been shown that we can successfully choose the correct model of rotational profile among the prepared models.
Then, we have applied the Bayesian scheme to one of the Kepler targets, KIC 11145123, for which the fast-convective-core rotation has been suggested by Hatta et al. (2019). Focusing on the convective-core rotation of the star, four models of rotational profile are constructed. The global likelihoods of the four models thus computed clearly indicate that the models with the fast-convective-core rotation is favored, supporting the previous suggestion by Hatta et al. (2019) from a Bayesian perspective. The estimated parameters have been used to construct an artificial rotational profile and the corresponding rotational shifts, based on which the validity of the obtained results have been checked.
In addition to the inference on the convective-core rotation, it has been suggested that the position of the rotational velocity shear is not at the convective boundary but located within the overshoot zone. Since it is generally considered that there are still numerous uncertainties in physics around the boundary between the convective core and the radiative region above for earlytype main-sequence stars (e.g., the position of the convective boundary, the extent of the overshoot zone, the possible rotational velocity shear, the possible dynamo mechanisms, and so on), the results of the study could pose a unique challenge to, for instance, numerical simulations of the dynamo mechanisms inside the convective core of early-type main-sequence stars.