© ESO 2022

1 Introduction

The past detection and characterization of exoplanets have provided rich information on their macroscopic properties such as orbital parameters, bulk properties (mass, radius, equilibrium temperature), and properties of their host stars. We are now moving from the discovery phase to a new era where we attempt to understand the atmospheric process and internal structure of exoplanets. Modeling the interior of exoplanets is currently one of the hot topics in exoplanetary sciences, yielding more findings than measurements of planetary mass and radius. Exoplanets are extremely diverse in the Universe, including terrestrial and oceanic Earth-sized planets, super-Earths, icy Neptune-like planets, and gas giants (Batalha 2014; Spiegel et al. 2014). The interior characterization of Earth-sized exoplanets mostly resorts to terrestrial bulk compositions from the knowledge of terrestrial planets in our Solar System (Sotin et al. 2007; Boujibar et al. 2020; Zhao & Ni 2021). Bulk refractory abundances from observations of host star’s photospheric compositions can be used to constrain the interior models of super-Earths and sub-Neptunes (Valencia et al. 2007; Dorn et al. 2017; Brugger et al. 2017; Adibekyan et al. 2021). For gas giants where direct measurements of their atmospheres are possible for transiting exoplanets, atmospheric observations show large diversities in effective temperature and most of these gas giants have no direct analogues in our Solar System in terms of physical and orbital parameters. Despite this, the internal structure and interior composition of gas giants have been explored within a similar model for Solar System giants (Nettelmann et al. 2010; Kramm et al. 2012). Thermal and structural evolution models have also been constructed to constrain the interior structures of gas giant exoplanets (Guillot & Showman 2002; Guillot 2005; Fortney et al. 2011; Thorngren et al. 2016; Fortney et al. 2021), which play a crucial role in understanding planet formation in diverse protoplanetary environments. Various interior models with different focuses have been developed for interior characterization. These models are capable of interpreting the measured mass and radius within their uncertainties. But they exhibit an inherent degeneracy issue as two planetary bodies with different internal structures or compositions could have the same mass and radius. In order to address the degeneracy issue, some inversion methods and new observations are proposed to constrain the interior structure of exoplanets. Rogers & Seager (2010) and Dorn et al. (2015, 2017) performed a full probabilistic inverse analysis of planetary compositions for rocky exoplanets with the uncertainties of measurements and models. Nettelmann et al. (2010) and Kramm et al. (2012) showed that the tidal Love number k₂ is a useful observable constraint on the interior models of gas giants, reducing the inherent degeneracy of interior models.

More recently, a wide range of machine learning techniques have been applied to Earth and planetary sciences, such as random forests (RFs), deep neural networks (DNNs), and mixture density networks (MDNs). The RF and DNN algorithms have powerful abilities in extracting information from high-dimensional data. In regression tasks, the objective of RFs and DNNs is to use some input variables to predict one or more variables, yielding discrete values for output variables. They have been widely applied to various fields in recent years, including exoplanet characterization (Ulmer-Moll et al. 2019), planetary growth models (Alibert & Venturini 2019), and the volcanic origin on Earth (Zhao et al. 2019). However, both RFs and DNNs show mediocre performance on degenerate problems where one set of input variables yields more than one solution for one output variable (Bishop 1994). The MDN algorithm is thus developed for degenerate cases without sacrificing degenerate solutions. It combines the structure of a conventional neural network with a probability mixture model, yielding multimodal probability distributions for output variables. MDNs have been shown by Atkins et al. (2016) to be a powerful method to infer the governing parameters for mantle convection. MDNs have also been shown by Baumeister et al. (2020) to be a useful tool to infer the distribution of possible thicknesses of each planetary layer for exoplanets with 0.01 M_⊕ < M < 25 M_⊕. Very recently, we trained MDNs to predict the radial structure and core properties of rocky exoplanets with 0.1 M_⊕ < M < 10 M_⊕ (Zhao & Ni 2021).This motivates further investigations of how the MDN works for understanding the interior of gaseous giant exoplanets. Gas giant exoplanets span over a much wider range of physical conditions than rocky exoplanets. In particular, they show larger diversities in thermal states. For example, their core properties and amounts of heavy elements are varied by several orders of magnitude and they are divided into cool, warm, and hot Jupiters according to thermal conditions. This brings about some difficulties for the trained MDN to accurately predict the interior properties of gas giant exoplanets. Nevertheless, once a MDN modelis well trained, users can quickly assess the interior structure of planets without calculating any additional interior models. In Sect. 2, we describe the interior model of gaseous giant exoplanets and present the practical aspects of internal structure calculations. Then we briefly introduce the MDN algorithm applied in this field in Sect. 3 since more details on MDNs have already been explained in the works of Baumeister et al. (2020) and Zhao & Ni (2021). Section 4 presents our results of the MDN training and testing, followed by the applications of the trained MDN to four representative gas giants. The predictive ability of the trained MDN is discussed in detail. Finally, the main results of this paper are summarized in Sect. 5.

2 Interior structure models and training data sets

The standard structure equations for a rotating body are as follows: $$\frac{\partial P}{\partial m}=-\frac{Gm}{4\pi r^4}+\frac{{\omega }^2}{6\pi r},$$(1) $$\frac{\partial r}{\partial m}=\frac{1}{4\pi r^2\rho },$$(2) $$\frac{\partial T}{\partial m}=\frac{T}{P}\frac{\partial P}{\partial m}{\nabla }_T,$$(3)

where ρ is the density at mean radius r, m is the mass enclosed within meanradius r, ω is the planet’s angular velocity, and∇_T is the temperaturegradient ∇_T ≡ (d lnT∕d lnP) which is equal to ∇_ad in the adiabatic case. The closure of the above system of Eqs. (1)–(3) requires the equation of state (EOS) of the H–He–rock–ice mixture, namely the density as a function of temperature, pressure, and composition. For hydrogen and helium, we chose to use the recent Chabrier et al. (2019) EOSs for pure hydrogen and helium, namely CMS19. The calculations of Chabrier et al. (2019) combined ab initio simulations for intermediate-density matter with semi-analytical calculations for low-density and high-density matter, spanning over a wide range of densities from 10⁻⁸ to 10⁶ g cm⁻³, with pressures from 10⁻⁹ to 10¹³ GPa and temperatures from 10² to 10⁸ K. The heavy elements, which can be divided into ices and rocks, are assumed to be either water for ices or serpentine for rocks (Thorngren et al. 2016). Here we adopt a 50–50 ice–rock mixture and use the Thompson (1990) ANEOS EOS for heavy elements. The EOS of the H–He–rock–ice mixture is generated using the separate EOSs for pure species in terms of the additive volume rule. The density for the mixture is therefore approximated at any given pressure P, temperature T, and composition (X, Y, and Z) as $$\frac{1}{\rho (P,T,Y,Z)}=\frac{X}{{\rho }^{\text{H}}(P,T)}+\frac{Y}{{\rho }^{\text{He}}(P,T)}+\frac{Z(1-f_{\text{rock}})}{{\rho }^{\text{ice}}(P,T)}+\frac{Z{f}_{\text{rock}}}{{\rho }^{\text{rock}}(P,T)},$$(4)

where the rock mass fraction in heavy elements f_rock is fixed at 0.50 as stated above.

For gas giants in our Solar System, the current knowledge from both observational and theoretical sides have demonstrated some inhomogeneities in their interiors, such as the hydrogen–helium immiscibility region in Saturn (Leconte & Chabrier 2013; Schöttler & Redmer 2018; Mankovich & Fortney 2020; Brygoo et al. 2021) and the dilute core region in Jupiter (Wahl et al. 2017; Guillot et al. 2018; Ni 2019). For a planet with an inhomogeneous layer inside, the inhomogeneous layer would reduce the outward heat transport, causing the deeper interior to get hotter and exhibit lower densities. Such a density deficit would be compensated by adding more heavy elements in the interior in order to comply with the measured mass and radius of the planet (Ni 2019; Nettelmann et al. 2021). Atmosphere grids, which serve as an upper boundary condition for interior models, have been computed through radiative–convective equilibrium atmosphere models to constrain the atmospheric structure and thermal evolution of gas giants (Fortney et al. 2007, 2011). Fortney et al. (2007) computed grids of radiative-convective model atmospheres for hydrogen–helium planets on a large grid of gravities, intrinsic effective temperatures, and incident fluxes. For a highly irradiated hot Jupiter or an old and relatively low-mass planet where incident stellar flux dominates over intrinsic flux, an atmospheric radiative region grows and quickly becomes almost isothermal (Guillot & Showman 2002; Fortney et al. 2007). The appearance of isothermal regions makes the deep interior adiabat cooler than the adiabatic case. This is indeed the opposite in regards to the appearance of inhomogeneous regions.

Atmosphere grids for Jupiter and Saturn were also computed by Fortney et al. (2011) over a range of intrinsic fluxes and surface gravities. These grids provide the temperature at the 10-bar level T₁₀ over a wide range of the intrinsic effective temperature T_int and surface gravity g. However, they could not cover all atmospheric conditions of giant planets along the evolutionary track. In order to cover more atmospheric conditions, Leconte & Chabrier (2013) made an extrapolation by an analytic fit T₁₀ (T_int, g) to the numerical grids of Fortney et al. (2011). Here, we adopted Jupiter’s analytic function T₁₀ (T_int, g) to infer the intrinsic temperature T_int of the modeled planet with a given temperature T₁₀ and gravity g. This procedure works well for Jupiter-like giants. For strongly irradiated planets such as hot Jupiters, Guillot & Showman (2002) showed the atmosphere grids suitable for isolated objects would overestimate the atmospheric temperature. In turn, the isolated case only provides a lower boundto T_int. Various atmospheric processes affecting the atmospheric temperature are very complex and beyond the scope of this work.

The uncertainties from internal structures and atmospheric grids are often overshadowed by large uncertainties in available observationsof mass, radius, and effective temperature for gas giant exoplanets. In this work, we do not attempt to consider all the modeling uncertainties, but grasp the key ingredients to make the generated data as realistic as possible. We modeled the internal structure of gaseous giant exoplanets using two-layer structure models. A gaseous exoplanet consists of an isothermal central compact core made of heavy elements and a homogenous convective envelope composed of hydrogen–helium mixtures with a certain amount of heavy elements. Following the interior models of gas giants in our Solar System, we define the surface of gaseous giant exoplanets as the 1-bar pressure level and set the temperature at the 1-bar level as T₁. The mass fractions of helium and heavy elements in the envelope are denoted by Y and Z, respectively. The Y∕(1 − Z) value was constrained so that the mean helium mass fraction of the planet would match the value in the protosolar nebula Y_proto∕(X_proto + Y_proto) = 0.277 ± 0.006 (Serenelli & Basu 2010; Guillot et al. 2018). The central compact core, which is full of heavy elements, has the mass of M_core and its thermal state is characterized by the pressure and temperature at the core-envelope boundary (CEB), P_CEB and T_CEB. Figure 1 illustrates a schematic diagram of the interior structure of gaseous giant exoplanets.

Giant planets in our Solar System serve as calibrators for the interior modeling of gas giant exoplanets. In generating the training data sets, we imposed the following constraints on our interior models based on our current knowledge of gas giants, especially the ones in our Solar System. We first focused our attention to gas giant exoplanets with masses M between 0.1 and 10 Jupiter masses (M_J), and the core mass M_core was randomly sampled between 0 and 50 M_⊕. Owing to planet rotations, a first-order correction to the gravitational potential appears in the hydrostatic equilibrium Eq. (1). Here, the rotation period of modeled exoplanets, hardly having an influence on planetary interiors, was set to a random value between 0.1 and 10.0 days. Many gas giant exoplanets are much closer to their host stars and hence much hotter than gas giants in our Solar System. So the surface temperature T₁ was randomly sampled between 100 to 3000 K. As for the heavy-element abundance Z in the envelope, it is known that Jupiter’s atmosphere exhibits the heavy-element enrichment that is about 3 times that of solar abundances and Saturn’s atmosphere is more enriched in heavy elements by a factor of roughly 10 relative to solar abundances. With this in mind, the heavy-element abundance Z is assumed to randomly vary between 1.0 and 20.0 times solar abundances. We note that these restrictions do not cover all gaseous giant exoplanets. After initializing the interior structure of modeled planets, we numerically integrated the system of Eqs. (1)–(3) until convergence. In total, we obtained ~200 000 data points for the training data set.

Fig. 1 Schematic interior structure of gaseous giant exoplanets. The compact core with mass Mcore is characterized by the pressure and temperature at the core envelope boundary (CEB), PCEB and TCEB. The surface defined by the 1-bar pressure level is characterized by the temperature T1 bar and the deeper atmosphere at the 10-bar pressure level by the temperature T10 bar.

Fig. 2 Sketch map of the mixture density network (MDN) where Gaussian mixture models are combined with the neural network. Three hidden layers are involved in the neutral network and each layer has 512 neurons. Twenty Gaussian mixture components are involved in the mixture model and each component is characterized by means μ, variances σ2, and mixing coefficients α.

3 Machine learning approaches

3.1 Mixture density networks

In most recent years, MDN technology, composed of neural networks and mixture models, has been successfully applied to characterize the interior structure of exoplanets (Baumeister et al. 2020; Zhao & Ni 2021), where the multimodal probability distribution instead of one discrete value is predicted for each output variable. Figure 2 illustrates the MDN structure with Gaussian mixture models combined with the neural network. Following Zhao & Ni (2021), three hidden layers with 512 neurons each were used. That is, a neural network of “Input–512–512–512–Output” was built. For the output layer, a linear association of 20 Gaussian mixture components was adopted for each output variable and each Gaussian mixture component was characterized by means μ, variances σ², and mixing coefficients α.

The data we generated using the above interior models are randomly divided into training and testing data sets by a ratio of nine to one. And 10% of the training data set was used as a validation data set to validate the model results at the end of each epoch. The Adaptive moment estimation (Adam) optimizer (Kingma & Ba 2014)with a learning rate of 0.001 was chosen in the training process and Rectified Linear Unit (ReLU) functions were used to activate the hidden layers. The mixture model layer uses logit functions for mixing coefficients α neurons and nonnegative exponential linear unit (NNELU) functions for μ and σ² neurons. Early-stopping (Montavon et al. 2012) and dropout (Srivastava et al. 2014) methods were employed in the MDN training to avoid the occurrence of overfitting. Min-max normalization was used in order to prevent introducing any bias in the model training. The difference between the MDN outputs and the target data was measured by the so-called loss function and the average negative log-likelihood (NLL) was adopted, $$\mathcal{L}(\Theta )=-\frac{1}{N}{\sum }_{(x,y)\in \mathbb{D}}\log \text{Bigg}[{\sum }_{i=1}^m{\alpha }_i(x) {\varphi }_i\text{Big}(y | {\mu }_i(x),{\sigma }_i^2(x)\text{Big})\text{Bigg}],$$(5)

where Θ represents the MDN model parameters, N is the size of the training data set $\mathbb{D}$ consisting of input vectors x and target vectors y, m is the number of the Gaussian mixture components, α_i is the mixture weight of the ith Gaussian distribution, ϕ_i represents the ith Gaussian distribution of the target vector y, and μ_i and ${\sigma }_i^2$ are the mean and variance of the ith Gaussian distribution. The iterative training process was performed by adjusting the weight values of the MDN until the NLL loss had a good convergence. Finally, the performance of the trained MDN model was assessed using the testing data set. More details about the MDN model can be found in Baumeister et al. (2020) and Zhao & Ni (2021).

Our MDN model was trained on NVIDIA GeForce GT 730 graphics card using Keras (Chollet et al. 2015) based on TensorFlow (Abadi et al. 2016) and MDN layers (Martin & Duhaime 2019). The scikit-learn library (Pedregosa et al. 2011) was used for data processing.

3.2 Model training and validation

Knowledge of the heavy-element mass of gas giants can provide important clues as to the mechanism of planet formation (Fortney et al. 2007; Thorngren et al. 2016). The heavy-element enrichment of Jupiter and Saturn, including their atmospheric compositions and core masses, have been extensively explored (Wahl et al. 2017; Guillot et al. 2018; Ni 2019; Militzer et al. 2019; Debras & Chabrier 2019; Mankovich & Fortney 2020; Ni 2020; Nettelmann et al. 2021), improving our understanding of planet formation. However, there is little knowledge of the bulk or atmosphericcompositions of exoplanets, providing much fewer constraints on planet formation from exoplanetary systems. In the MDN model, we pay particular attention to the heavy-element mass of gas giant exoplanets, $M_Z=M_Z^{\text{env}}+M_{\text{core}}$.

At present, direct measurements of an exoplanet’s atmosphere are possible for transiting gaseous exoplanets. One key diagnostic for a giant planet atmosphere is the equilibrium temperature T_eq, which represents the energy flux received from its host star. It is dependent on the luminosity of the star, the planet’s Bond albedo, and the distance between the planet and star. Meanwhile, the total effective temperature T_eff of a planet can be most easily derived from its experimentally measured infrared emission spectrum. The most unconstrained thermal quantity of a planet is the intrinsic effective temperature T_int, which varieswith the age and mass of a planet. We can express the total luminosity of a planet with radius R as $$L=4\pi R^2\sigma T_{\text{eff}}^4=4\pi R^2\sigma (T_{\text{eq}}^4+T_{\text{int}}^4).$$(6)

After remote measurements of the total emitted infrared flux and the reflected visible flux, the radius R could be constrained from the known T_int according to Eq. (6) (Marley et al. 2007). The intrinsic temperature T_int is also among the most important quantities in planetary evolution calculations (Guillot & Showman 2002; Fortney et al. 2007, 2011; Thorngren et al. 2016). So, the intrinsic temperature T_int is predicted as one MDN output as well. As for the thermodynamics properties of planetary cores, the pressure and temperature at the CEB, P_CEB and T_CEB, are predicted from the MDN similar to the case of rocky exoplanets (Zhao & Ni 2021).

We designed two comparative MDN models, with (M, R) inputs and with (M, R, T₁) inputs. First, we investigated how much data is required to train a good model by evaluating the performance of MDN fitting to the training data set of a different size. When the training data set reaches the size of 10⁴, as shown in Fig. 3a, the NLL loss changes little over the size for both MDN models. For the sake of comprehensiveness, we used all the training data to train our MDN models. Figure 3b and c illustrates the learning curves of the two MDNs, where the training and validation losses are shown as a function of the training epoch. The validation NLL loss reaches an approximate constant after 80 epochs for both MDN models. This implies that continuing training could not further improve the performance of the MDN model. Additionally, the MDN with (M, R, T₁) inputs yields the smaller NLL loss than the MDN with (M, R) inputs, suggesting its better training performance. After both MDN models were trained, the testing data set was fed to them to evaluate their performance in predicting the interior properties of gas giant exoplanets.

Fig. 3 Negative log-likelihood (NLL) loss for the MDN model trained on mass and radius inputs only and for the MDN model trained on mass, radius, and T1 inputs. (a) Dependence of the NLL loss on the size of the training data set for the two MDN models. (b) Learning curve of the MDN trained on mass and radius inputs only for both the training and validation data sets. (c) Learning curve of the MDN trained on mass, radius, and T1 inputs for both the training and validation data sets.

Fig. 4 Predicted distributions of four output variables from the MDN trained on (M, R) inputs versus the actual values obtained from the interior models for the testing data set: total heavy-element mass MZ, intrinsic effective temperature Tint, as well as pressure and temperature at the core envelope boundary (CEB), PCEB and TCEB. Blue dasheddiagonal lines denote a perfect performance of the MDN. The predicted distributions are colored according to the local probability density (the color scale is black at the maximum).

4 Results

Figure 4 illustrates the ability of the MDN trained on mass and radius inputs only to predict the interior properties of gaseous giant exoplanets, where the distributions predicted by the MDN from the testing dataset are displayed versus the actual values obtained from the interior models. With only mass M and radius R as inputs, the MDN can infer the CEB pressure P_CEB well. The similar behavior also emerged in the MDN predictions for rocky exoplanets (Zhao & Ni 2021). It is accounted for by a trade-off between the integral interval and pressure gradient of the hydrostatic Eq. (1). In contrast to the case of rocky exoplanets, the present MDN tends to overestimate the P_CEB values for Jupiter-sized planets and underestimate the P_CEB values for super-Jupiters. This can be easily understood because the present P_CEB values span over a much wider range from 10⁻² to 10⁴ Mbar, leading to some difficulties for the MDN to accurately predict the whole set of P_CEB. It is also becauseplanetary rotation affects the pressure profile of gas giants (see Eq. (1)), but the rotation period of a planet is not included as input to the MDN. It should be noted that these deviations from the actual P_CEB values can be systematically ignored with respect to the large uncertainties in available observations of planetary mass and radius. The CEB temperature T_CEB and intrinsic effective temperatureT_int are moderately constrained since the adiabatic temperature gradient is primarily a function of the hydrostatic pressure (see Eq. (3)). For the total mass of heavy elements, the MDN predictions exhibit a broad distribution, showing the considerabledegeneracy of possible interior structures.

The surface temperature T₁ is closely related to the thermal state of gas giants. It has been demonstrated by Baumeister et al. (2020) and Zhao & Ni (2021) that adding an additional input into the MDN could significantly improve the performance of the trained MDN. Therefore, the MDN was also trained on (M, R, T₁) inputs, compared with the MDN trained on mass and radius inputs only. Figure 5 illustrates the predictive ability of the MDN trained on (M, R, T₁) inputs. As can be seen, the total heavy-element mass M_Z was varied between 1 to 1000 M_⊕, the intrinsic temperature T_int between 50 and 2650 K, the pressure P_CEB between 10⁻² to 10⁴ Mbar, and the temperature T_CEB between 3 × 10³ to 5 × 10⁵ K. In view of such wide ranges, a small chart is inserted in order to illustrate the low-value behavior for each MDN output. Adding T₁ as an additional input significantly enhances the predictive ability of the MDN. The most considerable improvement lies in the heavy-element mass M_Z, implying that the heavy-element enrichment of a planet is closely correlated with planetary thermal states. This gives an active response to the findings for the interiors of Jupiter and Saturn (Ni 2019; Nettelmann et al. 2021). One can see that many giant exoplanets in the testing data set have the heavy-element mass of more than 50 M_⊕, suggesting significant amounts of heavy elements in the envelope rather than the compact core. This inference is only for the giant exoplanets considered in this work since the generated data set does not cover the overall sample of known giant exoplanets. The MDN also predicts T_int and T_CEB with a higher accuracy, where their distributions predicted by the MDN are almost concentrated along the diagonal. Using Jupiter’s interior as a reference, its four properties are marked in the inserted charts of Fig. 5. As one would expect, limited knowledge of H–He EOSs and phase separation creates some difficulties when calculating the interior structure of Jupiter and hence some uncertainties about its heavy-element mass and CEB temperature (Wahl et al. 2017; Guillot et al. 2018; Ni 2019; Debras & Chabrier 2019; Nettelmann et al. 2021). As shown from four inserted charts in Fig. 5, all four interior properties seem to be overestimated for small giants such as Jupiter-sized planets. When a planet is getting bigger, an underestimation of four interior properties appears instead and becomes more obvious. Such an unwelcome behavior can be attributed to a very wide range of four interior properties, as mentioned above.

Fig. 5 Same as in Fig. 4, but for the predicted distributions of four output variables from the MDN trained on (M, R, T1) inputs. In view of the wide ranges of these four variables, a small chart is inserted for each output variable to illustrate the low-value behavior. Jupiter’s interior properties are marked by navy blue stars in the four inserted charts.

Application to gas giants in the Solar System and gaseous giant exoplanets

We applied our trained MDN with (M, R, T₁) inputs to four representative planets: gas giants in our Solar System and beyond. Meanwhile, the validation data sets for these four planets were generated using the above interior models and the uncertainties on the mass and radius were considered by introducing a standard deviation $$\text{begin{align}}\sigma =\text{sqrt{big}[{(\Delta M/M)}^2+{(\Delta R/R)}^2\text{big}]\text{big}/2\text{}},\text{nn end{align}}$$(7)

where M and R are the actual values of planetary mass and radius, and ΔM and ΔR denote the deviations of one sample from the actual values.

Jupiter and Saturn, as gas giants in our Solar System, have been explored by various missions over the past 48 yr, providing strong constraints on their interior structures (Wahl et al. 2017; Guillot et al. 2018; Militzer et al. 2019). The basic thermal measurements of Jupiter and Saturn are given in Table 1. For direct comparison with the current interior model predictions (Wahl et al. 2017; Guillot et al. 2018; Ni 2019, 2020; Militzer et al. 2019; Nettelmann et al. 2021), the rotation periods of Jupiter and Saturn were separately fixed at 9 h 55 min 29.69s (Higgins et al. 1996) and 10 h 33 min 34 s (Militzer et al. 2019), and the standard deviation σ was restricted toa small value of 1% in the interior model calculations. The uncertainties on the surface temperatures T₁ was also taken into account in the interior model calculations, where T₁ was randomly varied between 160.0 and 170.0 K for Jupiter and between 130.0 and 140.0 K for Saturn. Figures 6 and 7 show the MDN predictions for Jupiter and Saturn using (M, R, T₁) inputs, compared with the validation data set consisting of the interior model solutions. In particular, we would like to note that both MDN predictions and interior model solutions are confined to a very narrow area with respect to the whole range of interior properties. For the sake of clarity, we plotted the distributions of the interior properties within a very narrow range, although this creates the illusion of large derivations of the MDN predictions from the interior model solutions. In order to examine whether theMDN is capable of describing the thermal effect on planetary interiors, we paid special attention to T₁ in the MDN predictions. Three different T₁ inputs of 160, 165, and 170 K were chosen for Jupiter with the mass and radius inputs unchanged and three different T₁ inputs of 130, 135, and 140 K for Saturn. The higher surface temperature T₁ seems to allow a larger amount of heavy elements in the planetary interior. This is quite consistent with the interior model results of Ni (2019). Also, as the surface temperature is increased, the MDN predictions move toward an area of larger T_int and hotter T_CEB.

Various interior models with new ingredients have been constructed to interpret the measured gravity data of Jupiter and Saturn. All interior models are dependent on the EOS of hydrogen and helium. Ni (2019) adopted different EOSs including CMS19 in four-layer structure calculations and investigated the effect of dilute cores on Jupiter’s interior. The CMS19 adiabat generally has a higher density (pressure) in the Mbar region with respect to the EOSs of Militzer & Hubbard (2013) and Saumon et al. (1995) (Chabrier et al. 2019). This would bring in a lower planetary metallicity for CMS19. Indeed, the results of CMS19 yield Jupiter’s total heavy-element mass over a narrow range of 8–12 M_⊕, showing the lowest global metallicity among four EOSs (Ni 2019). Militzer et al. (2019) explored Saturn’s interior from the Cassini Grand Finale gravity measurements using the MH13+SCvH EOS, where Saturn’s heavy-element mass is predicted as 16.5–23 M_⊕. Ni (2020) discussed the internal structure and interior composition of Saturn within the four-layer structure model using REOS3b and MH13+SCvH, where the MH13+SCvH solutions yield Saturn’s heavy-element mass of 16–19 M_⊕. With respect to the MH13+SCvH results, the models using CMS19 are expected to allow a smaller amount of heavy elements. That is, the CMS19 solutions would yield Saturn’s total heavy-element mass of roughly or less than 16 M_⊕. Very recently Nettelmann et al. (2021) have investigated the interiors of Jupiter and Saturn using the CMS19 H–He EOS, suggesting the heavy-element masses of 7.5–10.1 M_⊕ and 12.6–13.6 M_⊕ for Jupiter and Saturn. For comparison, the CMS19 results for the total heavy-element masses of Jupiter and Saturn, which were determined from the detailed interior models of Ni (2019) and Nettelmann et al. (2021), are also shown in Figs. 6 and 7. One can see that the MDN predictions of M_Z not only show good agreement with the validation data set within the error allowed, but also follow the other detailed interior models of Jupiter and Saturn (Nettelmann et al. 2021; Militzer et al. 2019; Ni 2019, 2020). There is an overestimation of T_int, P_CEB, and T_CEB for the MDN predictions, as discussed in the previous subsection. This overestimation can be neglected as compared with the systematic error because these quantities are varied by several orders of magnitude. Based on the thermal measurements listed in Table 1, the intrinsic temperatures of Jupiter and Saturn are derived from Eq. (6) as 98.2 K and 78.7 K, respectively. The MDN predictions of T_int cover these two values. The MDN predictions of P_CEB and T_CEB are also consistent with the four-layer interior models of Ni (2019, 2020).

For gaseous giant exoplanets, we consider CoRoT-9b as a Jupiter-sized representative and CoRoT-6b as a hot Jupiter representative. Observational data of these two exoplanets have been provided by several groups as summarizedbelow. Deeg et al. (2010) determined the absolute values for CoRoT-9b’s mass and radius as M = 0.84 ± 0.07M_J and R = 1.05 ± 0.04 Jupiter radii (R_J) and gave anirradiation temperature between 250 K and 430 K. Bonomo et al. (2017) performed a Bayesian combined analysis of photometric observations from the CoRoT and Spitzer space telescopes and five-year spectroscopic monitoring with HARPS. They obtained the best values for CoRoT-9b’s mass, radius, and equilibrium temperature as M = 0.84 ± 0.05 M_J, $R={1.066}_{-0.063}^{+0.075}$ R_J, and T_eq = 420 ± 16 K. Fridlund et al. (2010) analyzed the photometric light curve from the CoRoT telescope as well as the photometric and spectroscopic ground-based follow-up observations, resulting in the physical properties of CoRoT-6b as M = 2.96 ± 0.34M_J, R = 1.166 ± 0.035R_J, and T_eq = 1017 ± 19 K. Southworth (2011) also derived the physical properties of CoRoT-6b using homogeneous methods and by combining all available photometric data with measured spectroscopic parameters of host stars, M = 2.96 ± 0.39M_J, R = 1.185 ± 0.050R_J, and T_eq = 1025 ± 16 K, covering the previous findings of Fridlund et al. (2010). At present, there is little knowledge of surface temperatures for CoRoT-9b and CoRoT-6b. And, furthermore, the equilibrium temperature T_eq of a planet is neither an upper nor lower bound to its surface temperature (Guillot 2005). Considering the effects of the surface temperatures on interior properties had been illustrated in Figs. 6 and 7, we roughly assume the surface temperatures T₁ of CoRoT-9b and CoRoT-6b to be equal to their equilibrium temperature T_eq. Figures 8 and 9 show the MDN predictions for CoRoT-9b and CoRoT-6b using (M, R, T₁) inputs to the MDN, compared with the validation data set with the standard deviation σ = 3%. As one can see, the distributions of the interior properties are restricted into a narrow range. In both cases, the MDN predictions show good agreement with the interior model solutions, especially for CoRoT-6b.

Fig. 6 MDN predictions of the interior properties for Jupiter. Colored lines denote the combined Gaussian mixture prediction from the MDN with (M, R, T1) inputs, where different shades of color correspond to different surface temperatures (T1 = 160, 165, 170 K from light to dark). Shaded areas denote the validation results from the interior models with the standard deviation σ = 1%. As additional information, Jupiter’s total heavy-element mass, determined from the detailed interior models using the CMS19 H–He EOS as well (Ni 2019; Nettelmann et al. 2021), is shown by green vertical bands.

Fig. 7 Same as in Fig. 6, but for Saturn. For the MDN predictions denoted by colored lines, different shades of color correspond to different surface temperatures (T1 = 130, 135, 140 K from light to dark). Saturn’s total heavy-element mass, determined from the detailed interior models using the CMS19 H–He EOS as well (Nettelmann et al. 2021), is shown by green vertical bands.

Fig. 8 MDN predictions of the interior properties for CoRoT-9b. Colored lines denote the combined Gaussian mixture prediction from the MDN with (M, R, T1) inputs. Shaded areas denote the validation results from the interior models with the standard deviation σ = 3%.

Fig. 9 Same as in Fig. 8, but for CoRoT-6b.

5 Summary

As an extension of the previous work for rocky exoplanets (Zhao & Ni 2021), we applied MDN techniques to infer the interior properties of gaseous giant exoplanets with 0.1 M_J < M < 10 M_J. Compared with the case of rocky exoplanets, the crucial difficulty for gas giant exoplanets is that they show larger diversities in thermal states and their interior properties span over a much wider range by several orders of magnitude. We developed two-layer interior models using the H–He–rock–ice EOSs and generated sufficient data samples by solving the interior models with various physical parameters. The data were then trained based on two MDN models and the resulting MDNs are available on GitHub¹. The interior properties of a planet, as MDN output variables, contain the total heavy-element mass in planetary interior M_Z, the intrinsic effective temperatureT_int, and the pressure P_CEB and temperature T_CEB at the CEB. The MDN trained on mass and radius inputs only works in predicting the pressure P_CEB, but this yields the large degeneracy of the other interior properties in particular for M_Z, as shown in Fig. 4. Using the temperature at the 1-bar level T₁ as additional input, the new MDN shows good performance in predicting the interior properties of gas giant exoplanets, as shown in Fig. 5. This suggests that the surface temperature is necessary to characterize the interior of gas giant exoplanets owing to their large diversities in thermal states. In spite of the difficulties mentioned above, the MDN can be efficiently trained for gaseous giant exoplanets to bypass the interior model calculations. We also applied the well-trained MDN to our Solar System giants (Jupiter and Saturn) and two representative gas giants (CoRoT-9b and CoRoT-6b). For Jupiter and Saturn, the MDN predictions are not only in good agreement with the validation data set, but also consistent with the more sophisticated interiors models of Militzer et al. (2019), Ni (2019, 2020), and Nettelmann et al. (2021). The trained MDN is also found to work well for CoRoT-9b and CoRoT-6b, as shown in Figs. 8 and 9. The MDN offers a convenient and powerful tool for quick knowledge of the interiorsof gaseous giant exoplanets and hence more insight into planet formation in diverse protoplanetary environments.

Acknowledgements

This work is supported by the National Natural Science Foundation of China (Grant No. 12022517), the Science and Technology Development Fund, Macau SAR (File No. 0048/2020/A1 and 0002/2019/APD), and the Pre-Research Projects on Civil Aerospace Technologies of China National Space Administration (Grant No. D020308 and D020303).

References

All Tables

All Figures

	Fig. 1 Schematic interior structure of gaseous giant exoplanets. The compact core with mass Mcore is characterized by the pressure and temperature at the core envelope boundary (CEB), PCEB and TCEB. The surface defined by the 1-bar pressure level is characterized by the temperature T1 bar and the deeper atmosphere at the 10-bar pressure level by the temperature T10 bar.
In the text

	Fig. 2 Sketch map of the mixture density network (MDN) where Gaussian mixture models are combined with the neural network. Three hidden layers are involved in the neutral network and each layer has 512 neurons. Twenty Gaussian mixture components are involved in the mixture model and each component is characterized by means μ, variances σ2, and mixing coefficients α.
In the text

Fig. 3 Negative log-likelihood (NLL) loss for the MDN model trained on mass and radius inputs only and for the MDN model trained on mass, radius, and T1 inputs. (a) Dependence of the NLL loss on the size of the training data set for the two MDN models. (b) Learning curve of the MDN trained on mass and radius inputs only for both the training and validation data sets. (c) Learning curve of the MDN trained on mass, radius, and T1 inputs for both the training and validation data sets.

In the text

Fig. 4 Predicted distributions of four output variables from the MDN trained on (M, R) inputs versus the actual values obtained from the interior models for the testing data set: total heavy-element mass MZ, intrinsic effective temperature Tint, as well as pressure and temperature at the core envelope boundary (CEB), PCEB and TCEB. Blue dasheddiagonal lines denote a perfect performance of the MDN. The predicted distributions are colored according to the local probability density (the color scale is black at the maximum).

In the text

	Fig. 5 Same as in Fig. 4, but for the predicted distributions of four output variables from the MDN trained on (M, R, T1) inputs. In view of the wide ranges of these four variables, a small chart is inserted for each output variable to illustrate the low-value behavior. Jupiter’s interior properties are marked by navy blue stars in the four inserted charts.
In the text

Fig. 6 MDN predictions of the interior properties for Jupiter. Colored lines denote the combined Gaussian mixture prediction from the MDN with (M, R, T1) inputs, where different shades of color correspond to different surface temperatures (T1 = 160, 165, 170 K from light to dark). Shaded areas denote the validation results from the interior models with the standard deviation σ = 1%. As additional information, Jupiter’s total heavy-element mass, determined from the detailed interior models using the CMS19 H–He EOS as well (Ni 2019; Nettelmann et al. 2021), is shown by green vertical bands.

In the text

	Fig. 7 Same as in Fig. 6, but for Saturn. For the MDN predictions denoted by colored lines, different shades of color correspond to different surface temperatures (T1 = 130, 135, 140 K from light to dark). Saturn’s total heavy-element mass, determined from the detailed interior models using the CMS19 H–He EOS as well (Nettelmann et al. 2021), is shown by green vertical bands.
In the text

	Fig. 8 MDN predictions of the interior properties for CoRoT-9b. Colored lines denote the combined Gaussian mixture prediction from the MDN with (M, R, T1) inputs. Shaded areas denote the validation results from the interior models with the standard deviation σ = 3%.
In the text

	Fig. 9 Same as in Fig. 8, but for CoRoT-6b.
In the text