A Deep Learning Approach to Extracting Nuclear Matter Properties from Neutron Star Observations

Abstract

Understanding the equation of state of dense neutron-rich matter remains a major challenge in modern physics and astrophysics. Neutron star observations from electromagnetic and gravitational wave spectra provide critical insights into the behavior of dense neutron-rich matter. The next generation of telescopes and gravitational wave detectors will offer even more detailed neutron-star observations. Employing deep learning techniques to map neutron star mass and radius observations to the equation of state allows for its accurate and reliable determination. This work demonstrates the feasibility of using deep learning to extract the equation of state directly from observations of neutron stars, and to also obtain related nuclear matter properties such as the slope, curvature, and skewness of nuclear symmetry energy at saturation density. Most importantly, it shows that this deep learning approach is able to reconstruct realistic equations of state and deduce realistic nuclear matter properties. This highlights the potential of artificial neural networks in providing a reliable and efficient means to extract crucial information about the equation of state and related properties of dense neutron-rich matter in the era of multi-messenger astrophysics.

1. Introduction

The determination of the equation of state (EOS) of dense neutron-rich matter is one of the most significant and pressing challenges in modern physics and astrophysics [1,2,3]. Understanding the EOS has far-reaching implications for a broad range of phenomena, such as heavy-ion collisions, binary neutron star mergers, supernovae, and gravitational waves. Nuclear physics has been a primary field of research for studying the EOS, with numerous experiments and theoretical investigations carried out to probe the properties of dense nuclear matter [4,5,6,7,8,9,10,11,12,13]. Astrophysics also plays a crucial role in examining the EOS, with research efforts focusing on studying neutron stars and compact binary mergers involving at least one neutron star [14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35]. In recent years, numerous research facilities, such as telescopes and gravitational-wave detectors, have been established to enhance our knowledge and understanding of the EOS [1,2].

In cold neutron stars, the nucleonic part of the EOS can be written in terms of the energy per nucleon $E(\rho ,\delta )$, where $\rho$ is the total density, $\delta =({\rho }_{\mathrm{n}}-{\rho }_{\mathrm{p}})/\rho$ is the isospin asymmetry, and ${\rho }_{\mathrm{n}}$ and ${\rho }_{\mathrm{p}}$ represent the neutron and proton densities, respectively. In the parabolic approximation the EOS takes the form [36]

$$E(\rho ,\delta )=E_{SNM}\left(\rho \right)+E_{\mathrm{sym}}\left(\rho \right){\delta }^2,$$

(1)

where $E_{SNM}\left(\rho \right)$ is the energy per nucleon of symmetric nuclear matter (SNM) and $E_{\mathrm{sym}}\left(\rho \right)$ is the nuclear symmetry energy. However, the EOS of cold nuclear matter under extreme conditions is still a subject of intense debate and uncertainty, especially at supra-saturation densities. This uncertainty mainly arises due to the unknown high-density behavior of the nuclear symmetry energy $E_{\mathrm{sym}}\left(\rho \right)$ (see, e.g., Refs. [4,5] and references therein).

To extract the equation of state (EOS) of nuclear matter from first principles, one needs to solve quantum chromodynamics (QCD), the theory of strong interactions. However, model-independent results currently only span a narrow density regime. Specifically, for low densities near $1-2{\rho }_0$ (where ${\rho }_0=0.16$ fm${}^{-1}$ denotes the saturation density of symmetric nuclear matter), one can blend ab initio approaches with nuclear forces that stem from Chiral Effective Theory ($\chi$EFT), which yield well-defined error bounds [37,38,39,40,41,42,43,44,45,46]. For higher densities at $\rho \gtrsim 50{\rho }_0$, perturbative QCD calculations give reliable results [47,48,49,50,51,52,53]. However, at intermediate densities in the range of $2-10{\rho }_0$, reliable QCD predictions remain elusive [54]. Therefore, determining the EOS in this domain still relies on non-perturbative methods such as Monte Carlo simulation of QCD on a lattice (lattice QCD), which face formidable challenges in view of the sign problem in finite-density systems [55]. Accordingly, establishing the EOS at intermediate densities still hinges on phenomenological approaches that employ many-body techniques and effective interactions, such as density functionals based on Skyrme, Gogny, the relativistic mean field (RMF) theory, or Similarity Renormalization Group (SRG)-based methods.

Recent years have seen remarkable strides in determining the EOS at high densities, driven by both nuclear laboratory experiments and multi-messenger astrophysical (MMA) observations of neutron stars (NSs). Specifically, experimental data from heavy-ion reactions conducted at intermediate to relativistic energies, in particular various forms of nucleon collective flow and kaon production, have already imposed significant constraints on the EOS of symmetric nuclear matter up to approximately $4.5{\rho }_0$ [6]. Moreover, the nuclear and astrophysical communities have joined forces and have made impressive headway in constraining the symmetry energy around and below the saturation density of nuclear matter, leveraging a combination of terrestrial nuclear experiments and astrophysical observations [5,9,10,11,56,57,58,59]. Nevertheless, the characteristics of the nuclear symmetry energy $E_{sym}\left(\rho \right)$ at densities above saturation and the potential occurrence of a transition from hadronic matter to quark matter remain the aspects of the high-density equation of state that harbor the most significant uncertainties [4,5,21,22,24]. The high-density behavior of $E_{sym}\left(\rho \right)$ is also importnat for the occurrence of novel particles, such as hyperons and resonances [60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75].

The recent MMA observations of NSs have offered a unique opportunity to explore the high-density EOS. This represents an alternative way to independently extract the EOS by means of statistical approaches (as highlighted in Refs. [35,76,77,78,79,80,81]). The diverse set of methods employed in these investigations encompasses Shapiro delay measurements of massive ∼$2M_{\odot }$ pulsars [82,83,84], radius measurements of quiescent low-mass X-ray binaries and thermonuclear bursters [76,77,78,85,86], X-ray timing measurements from the Neutron Star Interior Composition Explorer (NICER) mission [87,88,89], and detection and inference of gravitational waves from compact binary mergers involving NSs [90,91,92] by the LIGO/VIRGO/KAGRA collaboration [93,94,95]. Neutron stars are characterized by common observables including mass M, radius R, moment of inertia I, quadrupole moment Q, dimensionless tidal deformability $\Lambda$ (and its derivatives such as Love number $k_2$ and tidal deformability $\lambda$), and compactness $M/R$. The NICER mission, for instance, is directed towards determining the compactness $M/R$ of neutron stars by observing the gravitational lensing effect of thermal emission from the star’s surface. On the other hand, the tidal disruption of the star in the presence of its companion is quantified through the tidal deformability parameter $\lambda$, which is inferred from gravitational-wave (GW) observations of binary neutron star (BNS) and neutron star-black hole (NSBH) mergers.

There exist various statistical approaches to determine the most likely EOS from neutron star observational data. Of these, the use of Bayesian inference is widespread [35,76,77,78,79,80]. Gaussian processes also provide a non-parametric representation of the EOS [81]. However, the uncertainty in Bayesian analyses raises questions regarding the true nature of the dense matter EOS [54]. The major sources of uncertainties in Bayesian inference methods are due to the choice of prior probabilities, which represent subjective beliefs or assumptions about the parameters before incorporating data. Different choices of priors can lead to different results and interpretations [54]. In light of this, alternative model-independent methods are being sought. Deep neural networks (DNNs) [96,97] have garnered attention in the research community, where deep learning (DL) algorithms have displayed exceptional proficiency in tasks such as image recognition [98] and natural language processing [99]. Furthermore, these techniques have been applied in various physics and astrophysics domains, including the analysis of GW data for detection [100,101,102,103,104,105,106,107], parameter estimation [108,109], and denoising [110]. In previous works, we employed Convolutional Neural Network (CNN) [111] algorithms to detect and infer GW signals from BNS [112,113] and, very recently, from NSBH [114] mergers. Additionally, the use of DNNs as a tool to extract the dense matter EOS from neutron star observations has also been explored in a growing number of studies [54,115,116,117,118,119,120].

DL methods offer important advantages over Bayesian approaches in several important aspects. One key advantage is the ability of DL algorithms to automatically learn complex patterns and features from large datasets without requiring explicit assumptions or prior knowledge. This allows DL algorithms to capture and extract complex nonlinearities and correlations that may not be easily discernible using conventional Bayesian techniques. It is precisely this unique intrinsic feature of the DL approaches that opens an alternative way for a model-independent determination of the dense matter EOS directly from neutron-star astronomical data.

In a recent investigation [121], we presented an innovative approach to determine the nuclear symmetry energy, $E_{sym}\left(\rho \right)$, by utilizing DL techniques in conjunction with astronomical observations of neutron stars. Our results demonstrate that deep neural networks have the capacity to accurately extract $E_{sym}\left(\rho \right)$ from a set of $M-R$ or $M-\Lambda$ NS observations. This approach offers a promising avenue for exploring the high-density behavior of $E_{sym}\left(\rho \right)$, which remains a challenging task in nuclear physics.

Motivated by the results of our previous work [121], in this paper, we further develop our DL approach and apply it to the problem of determining the EOS of dense neutron-rich matter and related nuclear matter properties from mass-radius $M\left(R\right)$ measurements of neutron stars. In particular, we pay special attention to deducing the slope, curvature, and skewness of the nuclear symmetry energy, in addition to the EOS. Our results demonstrate that DL algorithms can accurately and reliably extract the NS EOS and nuclear matter properties from observational data. Most importantly, we show that our DL framework is able to successfully reconstruct realistic EOSs and nuclear matter properties, which brings us one step closer to revealing the true nature of the EOS of dense neutron-rich matter.

In the present work, we have structured our discussion as follows. In the Section 1, we have provided a brief introduction. In Section 2, we have presented the main aspects of our formalism. This includes a comprehensive overview of the essential features of the EOS employed in our analysis, along with the details of our DL algorithms, such as data generation, neural network architectures, and training methodologies. Subsequently, in Section 3, we have presented our results and their implications. Finally, in Section 4, we have summarized our findings and provided future research directions.

Conventions: We use units in which $G=c=1$.

2. Formalism

In this section, we introduce the methods employed in our study. First, we give a comprehensive overview of the main characteristics of the EOS used in the present investigation. Afterwards, we briefly recall the procedure for solving the structure equations of static neutron stars. In Section 2.3, we discuss the DL methodology employed in mapping the mass-radius $M\left(R\right)$ observations of neutron stars to the EOS, and also the procedure of mapping the reconstructed EOS to selected nuclear matter properties.

2.1. Equation of State

The determination of NS properties such as mass M and radius R heavily relies on the equation of state. In general, two theoretical approaches are used to determine the nuclear matter EOS—phenomenological and microscopic. Phenomenological approaches employ effective interactions that describe the ground state of finite nuclei. Approaches based on Skyrme interactions [122,123] and relativistic mean-field (RMF) models [124] are examples of such methods that have been extensively used in the study of low-density nuclear systems. However, their suitability for describing high isospin asymmetry systems and high-density regimes is limited, and without experimental data to constrain these models under such circumstances, their predictions may not be reliable [125]. On the other hand, microscopic approaches use realistic two-body and three-body interactions to describe the behavior of nucleons. Such realistic nuclear interactions are based on meson-exchange theory [126,127] or more recent $\chi$EFT [42,128,129,130]. Examples of microscopic many-body methods include the Brueckner–Hartree–Fock (BHF) approach [131], the Dirac–Brueckner–Hartree–Fock (DBHF) approach [132,133], the variational method [134], the Quantum Monte Carlo technique and its derivatives [135,136], the self-consistent Green’s function technique [137], $\chi$EFT [45], and the $V_{low;k}$ approach [138]. The primary difficulty encountered by these methods is how to handle the repulsive core of the nucleon–nucleon interaction at short ranges, and the way they approach this differentiate them from each other.

The nucleonic component of the EOS can be described by two quantities: the binding energy of symmetric nuclear matter, $E_{SNM}\left(\rho \right)$, and the symmetry energy, $E_{sym}\left(\rho \right)$ (see Equation (1)). These two quantities can be expanded as Taylor series around ${\rho }_0$ as:

$$E_{SNM}\left(\rho \right)=E_0+\frac{K_0}{2}{x}^2+\frac{J_0}{6}{x}^3,$$

(2)

$$E_{sym}\left(\rho \right)=S_0+Lx+\frac{K_{sym}}{2}{x}^2+\frac{J_{sym}}{6}{x}^3,$$

(3)

where $x\equiv (\rho -{\rho }_0)/3{\rho }_0$. The coefficients of these expansions can be related to various physical properties of nuclear matter and can be experimentally constrained. They have the following meanings [139]: $E_0\equiv E_{SNM}\left({\rho }_0\right)$, $K_0\equiv {[9{\rho }^2{d}^2{E}_{SNM}/d{\rho }^2]}_{{\rho }_0}$, and $J_0\equiv {[27{\rho }^3{d}^3{E}_{SNM}/d{\rho }^3]}_{{\rho }_0}$ are the binding energy, incompressibility, and skewness of SNM; $S_0\equiv E_{sym}\left({\rho }_0\right)$, $L\equiv {[3\rho d{E}_{sym}/d\rho ]}_{{\rho }_0}$, $K_{sym}\equiv {[9{\rho }^2{d}^2{E}_{sym}/d{\rho }^2]}_{{\rho }_0}$, and $J_{sym}\equiv {[27{\rho }^3{d}^3{E}_{sym}/d{\rho }^3]}_{{\rho }_0}$ are the magnitude, slope, curvature, and skewness of the nuclear symmetry energy at saturation density. Currently, the most likely values of these coefficients are known within certain ranges: $E_0=-15.9\pm 0.4$ MeV, $K_0=240\pm 20$ MeV, $-300\le J_0\le 400$ MeV, $S_0=31.7\pm 3.2$ MeV, $L=58.7\pm 28.1$ MeV, $-400\le K_{sym}\le 100$ MeV, and $-200\le J_{sym}\le 800$ MeV; as reported in Ref. [140]. As seen above, several of these parameters have rather moderate uncertainty. For example, the binding energy $E_0$ is estimated to be $-15.9\pm 0.4$ MeV, while the magnitude of the symmetry energy $S_0$ is $31.7\pm 3.2$ MeV. However, many of them still have significant uncertainty, such as the curvature of the symmetry energy $K_{sym}$, which could range from $-400$ MeV to 100 MeV, and the higher order coefficients, $J_0$ and $J_{sym}$, with even wider uncertainty ranges.

Although the Taylor expansions described by Equations (2) and (3) are known to diverge at higher densities [141], they can also be interpreted as parameterizations with free parameters [140]. This dual nature implies that, for systems with low isospin asymmetries, the Taylor expansions are valid near saturation density, while for highly neutron-rich systems at supra-saturation densities, Equations (2) and (3) should be treated as parameterizations [140]. More details on the relationship between the Taylor expansions and the parameterizations can be found in Ref. [140]. These expressions are frequently used to model the EOS of neutron stars and have been used, for instance, to solve the inverse structure problem of neutron stars, to constrain the high-density symmetry energy through neutron star observational data [140,142], and to determine the most likely values of high-density EOS parameters through Bayesian analyses using neutron star data [143]. Compared to the widely used piecewise polytropes, these parameterizations have the advantage of including isospin dependence and composition information throughout the density range, while still allowing for modeling a wide range of EOSs as predicted by various many-body approaches. This characteristic is particularly important for deducing the high-density behavior of $E_{sym}\left(\rho \right)$, as the parameterizations clearly separate the contribution of the symmetry energy to the EOS. These parameterizations were also instrumental in our previous work [121], in which we demonstrated the feasibility of deep neural networks to extract the nuclear symmetry energy directly from neutron-star observational data.

In this work, we employ the EOS metamodel described above to facilitate the procedure of extracting the EOS and selected nuclear matter properties from NS observational data. By varying the EOS parameters, we generated a multitude of EOSs and the associated $M-R$ sequences through the solution of the NS structure equations. The matter in the core of the neutron star is modeled as a mixture of neutrons, protons, electrons, and muons in beta-equilibrium (referred to as the $npe\mu$-model). We use Equations (2) and (3) to calculate $E_{SNM}\left(\rho \right)$ and $E_{sym}\left(\rho \right)$, and then evaluate $E(\rho ,\delta )$ through Equation (1). The pressure of the neutron star matter in $\beta$-equilibrium is obtained from the energy density, $\epsilon (\rho ,\delta )=\rho [E_n(\rho ,\delta )+M_N]+{\epsilon }_l(\rho ,\delta )$, where $M_N$ is the mean nucleon mass and ${\epsilon }_l(\rho ,\delta )$ is the lepton energy density, as

$$P(\rho ,\delta )={\rho }^2\frac{dϵ(\rho ,\delta )/\rho }{d\rho }.$$

(4)

The computation of ${\epsilon }_l(\rho ,\delta )$ is explained in detail in Ref. [144]. When the density of the neutron star matter falls below approximately 0.07 fm${}^{-3}$, the core EOS is complemented by a crustal EOS that is more appropriate for lower density regions. Specifically, we adopt the EOS by Pethick et al. [145] for the inner crust and the EOS by Haensel and Pichon [146] for the outer crust.

In this analysis, we use Equations (2) and (3) as parameterizations together with the parabolic approximation of the nucleonic EOS, given by Equation (1). We keep $E_0$ and $S_0$ constant and vary the rest of the parameters. The fixed values of $E_0$ and $S_0$ were obtained through a combination of theoretical calculations and nuclear laboratory experiments, which represent their most probable current values. By varying the remaining parameters, $K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$, we were able to generate multiple samples of the EOS. Our expanded parameter space compared to our previous work [121] allowed us to model a wider range of EOSs predicted by various many-body approaches and models of the nuclear interaction. The effect of varying these parameters is shown in Figure 1. Although these parameters are, in principle, free, their ranges are constrained by the asymptotic boundary conditions of the EOS near ${\rho }_0$ and $\delta =0$, as well as the requirement that the EOSs satisfy causality and the microscopic stability condition. Furthermore, the resulting NS models must be capable of supporting a maximal mass of at least 2.14M${}_{\odot }$, which is the heaviest observed pulsar so far [84]. Figure 2 displays the ranges of $E_{SNM}$ and $E_{sym}\left(\rho \right)$ that satisfy all of the constraints.

2.2. Structure Equations of Static Neutron Stars

In this section, we briefly revisit the procedure for calculating the mass M and radius R of static neutron stars. In the case of a spherically symmetric relativistic star, the Einstein field equations can be simplified to the Tolman–Oppenheimer–Volkoff (TOV) equation [147], as follows:

$$\frac{dP\left(r\right)}{dr}=-\frac{\epsilon \left(r\right)m\left(r\right)}{r^2}\left[1+\frac{P\left(r\right)}{\epsilon \left(r\right)}\right]\left[1+\frac{4\pi r^{3}P\left(r\right)}{m\left(r\right)}\right]{\left[1-\frac{2m\left(r\right)}{r}\right]}^{-1},$$

(5)

where the mass within a sphere of radius r is determined by

$$\frac{dm\left(r\right)}{dr}=4\pi \epsilon \left(r\right)r^2.$$

(6)

To solve the above equations, one needs to supplement them with the EOS in the form $P\left(\epsilon \right)$. The integration of Equations (5) and (6) starts at the center of the star ($r=0$) with the initial conditions $m(r=0)=0$ and ${\epsilon }_c=\epsilon (r=0)$. Some care should be taken at $r=0$ since the above equations are singular at the center. The integration continues until the pressure P reaches zero, indicating the boundary of the star at $r=R$. The radius and gravitational mass of the neutron star are then given by R, where P vanishes, and $M=m\left(R\right)=4\pi {\int }_0^R\epsilon \left(r^{\prime }\right)r^{\prime 2}d{r}^{\prime }$, respectively.

For a given EOS, there is a unique relationship between the stellar mass and the central density ${\epsilon }_c$. Thus, for a particular EOS, there is a unique sequence of NSs parameterized by the central density (or equivalently the central pressure $P_c=P\left(0\right)$). In Figure 3 we show the range of possible EOSs ($P\left(\rho \right)$) satisfying all constraints (left window), and the resultant $M-R$ NS sequences (right window).

2.3. Artificial Neural Networks

This section provides a brief overview of the setup, structure, and workflow involved in implementing deep neural networks (DNNs) for our particular application. Readers interested in a more in-depth discussion are encouraged to consult relevant machine learning articles [96,148] and textbooks [97,149].

We used two DNNs that have similar architectures to accomplish two different tasks: extracting the EOS of dense neutron-rich matter from a set of mass-radius neutron-star measurements and inferring specific nuclear matter properties from the resulting NS EOS in $\beta$-equilibrium. We refer to the first neural network EOS DNN (equation of state deep neural network) and the second one NuPRO DNN (nuclear properties deep neural network). Figure 4 provides a schematic representation of the procedure of using these DNNs to extract the EOS and selected nuclear matter properties.

Before we discuss the details of the deep learning methods employed in this work, we reiterate briefly the key benefits of using DL techniques for extracting the EOS of dense matter and nuclear matter properties directly from neutron star observational data. The major advantage of using DL is the ability of the trained models to extract more accurate and unbiased information about the EOS from the NS measurements. In contrast to conventional approaches that rely on pre-existing models and assumptions regarding the underlying physics, deep learning techniques can directly learn the complex relationships between the NS observables and the EOS. This data-driven paradigm allows for a model-independent extraction of the EOS, where the DL algorithms are able to capture highly nonlinear and intricate features in the data, leading to accurate and unbiased determinations of the EOS. Furthermore, compared to conventional methods, such as Bayesian inference, deep learning approaches can be significantly faster, making them a highly efficient tool for exploring and extracting key insights from large data sets.

2.3.1. EOS Network (EOS DNN)

In this study, we employed a supervised DL method to obtain the EOS by formulating a regression problem. The DNN takes as input $M\left(R\right)$ sequences, which are sets of points representing pairs of NS mass-radius measurements, and produces estimates of the EOS ($P\left(\epsilon \right)$) as output. Accordingly, the datasets used for the training, validation, and testing of EOS DNN consist of $M\left(R\right)$ sequences and $P\left(\epsilon \right)$ samples. We apply the EOS model introduced in Section 2.1 and vary the parameters in Equations (2) and (3) to generate many samples of the EOS, and subsequently by solving the NS structure equations, the corresponding $M-R$ sequences. Specifically, we set $E_0=15.9$ MeV and $S_0=31.7$ MeV, and vary the rest of the parameters by randomly sampling their values from their respective ranges: $K_0=240\pm 20$ MeV, $-300\le J_0\le 400$ MeV, $L=58.7\pm 28.1$ MeV, $-400\le K_{sym}\le 100$ MeV, $-200\le J_{sym}\le 800$ MeV. The latest results from the PREX collaboration have suggested a high upper limit of 143 MeV for the value of L [150]. However, we have not examined the effects of higher L values in this study and leave it to future work. We ensured that the resulting EOSs $P\left(\epsilon \right)$ satisfy two constraints: (i) the microscopic stability condition, $\frac{dP}{d\epsilon }\ge 0$, and (ii) the causality condition, $c_s\equiv \sqrt{\frac{dP}{d\epsilon }}\ge c$, where c is the speed of sound. In addition, the resultant NS models must be able to sustain a maximal mass of at least 2.14 M⊙[84]. These constraints limit the possible values of $K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$ and result in a final set of EOS samples. To generate simulated observational data for NS, we randomly selected 50 points from a genuine $M-R$ sequence, within the range of 1M${}_{\odot }$ to M${}_{max}$ allowed by the given EOS. Then, each input sample is an array of dimension $2\times 50$ consisting of 50 pairs of (M, R) values. The values of M and R are scaled by dividing them by 3 and 20 respectively to ensure that the input data are in the $(0,1)$ range. Similarly, each output sample is an array of dimension $2\times 50$ consisting of 50 pairs of estimated (P, $\epsilon$) values, representing the EOS in the density range from ∼0.4 ${\rho }_0$ to 5 ${\rho }_0$. The DNN in this context is responsible for mapping an input sequence, which represents pairs of NS mass-radius measurements, to an output sequence that represents the EOS, $P\left(\epsilon \right)$.

Here we need to mention that as of now, the prospect of obtaining a significant number of simultaneous mass and radius measurements of neutron stars might seem overly optimistic. Nevertheless, the rapid advancements in the development of next-generation telescopes and gravitational wave detectors hold the promise of a substantially larger number of neutron star observations in the near future.

In supervised learning, the data are divided into training, validation, and testing data sets. The training data set is used by the DNN to learn from, the validation data are used to verify whether the network is learning correctly, and the testing data are used to assess the performance of the trained model. Here, the training dataset consist of 120,000 independent $M\left(R\right)$ sequences, representing the DNN inputs, and 120,00 matching EOS samples, $P\left(\epsilon \right)$, representing the DNN outputs. From each $M\left(R\right)$ sequence we further draw 50 ensembles, each containing 50 randomly selected (M, R) pairs. In this way, each EOS sample in the training data set is represented by 50 different random ensembles drawn from the same genuine $M\left(R\right)$ curve. The final training data set therefore consist of $6\times {10}^6$ samples. Similarly, the final validation data set consist 250,000 samples, where each of 5,000 independent output EOS samples is represented by 50 different random ensembles drawn from the same $M\left(R\right)$ sequence. Finally, the testing data set consist of 5000 unique input and output samples, not used in the training and validation of the DNN.

EOS DNN is a feedforward neural network with 5 hidden, dense, fully connected layers of dimension 500, and ReLU activation functions, followed by a dense linear layer of dimension 100. The first layer corresponds to the input to the neural network, which, in this case, is a $2\times 50$ array containing the NS M and R values drawn from a given $M\left(R\right)$ curve. At the end, there is a linear output layer of dimension $2\times 50$ returning the estimated EOS, $P\left(\epsilon \right)$ (pairs of (P, $\epsilon$) points). The design of the network was optimized by fine-tuning several hyper-parameters such as the number and type of layers, the number of neurons in each layer, and the type of activation function. The optimal architecture was determined after multiple experiments and hyper-parameter tuning. A summary of the EOS DNN architecture is given in

Table 1. The EOS DNN (equation of state deep neural network) architecture comprises an input layer followed by 6 dense fully connected layers. The output layer returns the estimated equation of state, $P\left(\epsilon \right)$. The model has 1,102,600 trainable parameters. Further details can be found in the text.

	Layer	Activation	Size
	Input	–	2 × 50
1	Flatten	–	100
2	Dense	ReLU	500
3	Dense	ReLU	500
4	Dense	ReLU	500
5	Dense	ReLU	500
6	Dense	ReLU	500
7	Dense	Linear	100
8	Reshape	–	2 × 50
	Output	–	2 × 50

.

We used the Keras Python toolkit (https://keras.io (accessed on 25 March 2023)) which provides an easy-to-use API for the TensorFlow [151] (https://www.tensorflow.org (accessed on 25 March 2023)) deep learning library to build and train our neural network. The training was performed using stochastic gradient descent with adaptive learning rate with the ADAM method [152] and the AMSgrad modification [153]. For training, we used a batch size of 1000 and an initial learning rate of 0.001. We set a maximum of 5000 epochs for each training session, or until the validation error stopped decreasing. The best performing model was chosen using checkpoint-callback by selecting the one with the lowest loss value on the validation dataset. The training was executed on an NVIDIA Tesla V100 GPU, and the mini-batch size was determined automatically based on the characteristics of the GPU and data sets.

The cost or loss function used in our analysis was the mean absolute error (MAE), which is calculated as the average of the absolute differences between the predicted output values of the DNN model (${\widehat{y}}_i$) and the “true” output values ($y_i$) for each sample. The equation for the MAE is given by:

$$MAE=\frac{1}{n}\sum _{i=1}^n|{\widehat{y}}_i-y_i|,$$

(7)

where n represents the number of samples.

2.3.2. Nuclear Matter Properties Network (NuPRO DNN)

Once the $\beta$-stable matter EOS becomes available, we can proceed with the extraction of selected nuclear matter properties. In order to achieve this goal, we have trained another DNN, which we refer to as the NuPRO DNN. The input to the DNN consists of EOS samples represented as $P\left(\rho \right)$, which are sets of 50 equally spaced points within the interval $\rho =[0.08-0.8]f{m}^{-3}$, where the input data was converted to decimal logarithm values. On the other hand, the output corresponds to estimated selected nuclear matter properties, with respect to each input EOS sample. Our aim is to learn the mapping $\mathbf{y}\left(\mathbf{x}\right)$ with ${\mathbf{x}}_{\mathbf{i}}=[P_{\beta }\left({\rho }_1\right),P_{\beta }\left({\rho }_2\right),\dots ,P_{\beta }\left({\rho }_{50}\right)]$ representing the EOS and ${\mathbf{y}}_{\mathbf{i}}=[K_0,J_0,L,K_{sym},J_{sym}]$ being the corresponding set of parameters.

In this work, we used a training dataset comprising 120,000 samples of the equation of state $P\left(\rho \right)$, each with an accompanying set of parameters ($K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$) matching the specific EOS realization. Additionally, we constructed separate validation and testing datasets, each consisting of 5,000 data samples. The architecture of the NuPRO DNN model that we have implemented involves a feedforward structure with five dense hidden layers with dimensions 200, 200, 200, 100, and 50, respectively. We have used the ReLU activation functions for each of these layers. The neural network architecture of the EOS DNN was optimized through an iterative process involving multiple experiments and tuning of the hyper-parameters. The neural network’s input layer has a dimension of 50 and corresponds to the 50 uniformly distributed data points representing the equation of state, $P\left(\rho \right)$, while the output layer has a dimension of 5, which returns the estimated nuclear matter parameters. A summary of the network architecture is provided in

Table 2. The NuPRO DNN (nuclear matter properties deep neural network) architecture comprises an input layer with 50 dimensions, corresponding to the 50 equally spaced points of the equation of state $P\left(\rho \right)$. This is followed by five dense fully connected layers of varying dimensions, culminating in an output layer that returns the estimated nuclear matter parameters $K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$. The total number of trainable parameters in this deep neural network model is 118,555. Further information about this architecture can be found in the text.

	Layer	Activation	Size
	Input	–	50
1	Dense	ReLU	200
2	Dense	ReLU	200
3	Dense	ReLU	200
4	Dense	ReLU	100
5	Dense	ReLU	50
	Output	–	5

.

We used Keras and TensorFlow to develop and train our neural network. As before, we applied the ADAM method [152] with the AMSgrad modification [153] to implement stochastic gradient descent with an adaptive learning rate. In training the DNN, we chose a batch size of 1000 and set the initial learning rate to 0.001. We also restricted the number of epochs per training session to a maximum of 5000 or until the validation error stopped decreasing. We employed a checkpoint-callback to identify the model that achieved the lowest loss value on the validation dataset. The objective function chosen for this task was the mean-squared error (MSE), which measures the average of the squared differences between the predicted values of the DNN model, denoted by ${\widehat{y}}_i$, and the actual or "true” values, represented by $y_i$, over all the available samples, n. More specifically, the MSE is calculated as the sum of the squared differences divided by the number of samples, n, as shown in the below equation:

$$MSE=\frac{1}{n}\sum _{i=1}^n{({\widehat{y}}_i-y_i)}^2.$$

(8)

As with any numerical approach, DL techniques are not immune to uncertainties that can impact the reliability of their results. One significant source of uncertainty is the quality and quantity of data used to train the DNN. Insufficient or biased data may lead to erroneous or biased predictions. In addition, uncertainties could result from the different DNN architectures and hyperparameters used, such as the type and number of layers, activation functions, and learning rates. Different choices of neural network architecture and hyperparameters could cause different model outcomes. Moreover, model uncertainties could arise from the generalization ability of deep learning algorithms, as they may be affected by overfitting or underfitting, especially when presented with small datasets. Furthermore, DL methods can be sensitive to variations in input data, further contributing to uncertainties in their predictions. Therefore, it is critical to carefully consider and address these uncertainties to ensure the reliability and accuracy of DL predictions. We investigated the effect of “noise” in the DNN input data and the corresponding uncertainties in the DNN predictions in the next section.

3. Results

3.1. Extracting the EOS

We start by evaluating the ability of the DNN to reconstruct the EOS, $P\left(\epsilon \right)$, using a set of $M-R$ measurements that could result from electromagnetic observations of neutron stars, such as those from the NICER mission, for example. To accomplish this, we employ the previously described EOS DNN to a test dataset, which comprises ∼5000 simulated $M\left(R\right)$ sequences, and compare the resulting estimated output EOS with the exact EOS for each sample. Figure 5 displays the results for five representative examples from the test dataset, showing that the EOS (broken colored lines) for each input $M-R$ sequence agrees almost precisely with the "true” EOS (solid black lines) over the entire density range considered in this study. For visualization purposes, the EOS is presented in the form $P\left(\rho \right)$. The results for the remaining test data samples are comparable to those shown in Figure 5. Quantitatively, the mean absolute error across the entire test dataset is 0.5 MeV fm${}^{-3}$, with a standard deviation of 1.3 MeV fm${}^{-3}$. Choosing different sets of randomly-selected points from the authentic $M\left(R\right)$ curves has little effect on the accuracy with which the EOS is estimated.

A natural consequence of realistic observations of neutron stars is that they come with inherent uncertainties, which propagate into uncertainties in the estimated equation of state. To investigate the effect of the observational uncertainties on extracting the EOS via our DL approach, we prepared a test dataset assuming that the NS mass and radius measurements are subject to a measurement error. After randomly selecting N points ($M_i$, $R_i$) from a genuine $M\left(R\right)$ curve, with $\{i=1,2,\dots ,N;N=50\}$, we draw the actual M and R values from normal distributions, $\mathcal{N}({\mu }_M,{\sigma }_M)$ and $\mathcal{N}({\mu }_R,{\sigma }_R)$, respectively. In particular, we investigated the impact of smaller and larger observational errors on the estimated EOS. To simulate smaller uncertainties, we selected ${\sigma }_M=0.02$ M${}_{\odot }$ and ${\sigma }_R=0.1$ km, where ${\mu }_M=M_i$ and ${\mu }_R=R_i$, while to model larger uncertainties, we used ${\sigma }_M=0.1$ M${}_{\odot }$ and ${\sigma }_R=1$ km. To quantify the effect of observational uncertainties for each sample in the test dataset, for each case we drew 100 ensembles from the respective normal distribution and calculated the mean absolute errors. In Figure 6, we illustrate the effect of measurement errors for a representative example from the test dataset. It is seen that for smaller observational uncertainties (${\sigma }_M=0.02$ M${}_{\odot }$, ${\sigma }_R=0.1$ km) the estimated EOS $P\left(\rho \right)$ (broken green line) matches almost exactly the “true” EOS (solid black line). The greenish shaded band represents the corresponding mean absolute errors. For larger observational uncertainties (${\sigma }_M=0.1$ M${}_{\odot }$, ${\sigma }_R=1$ km), we see that at higher densities, of $\sim \rho /{\rho }_0\ge 4$, the reconstructed EOS (broken red line) starts to moderately diverge from the ground-truth EOS, but it is still within the reconstruction errors represented by the reddish shaded band (corresponding to the mean absolute errors). The trend observed in the presented results holds for the remaining test data samples as well, indicating that the precision of the EOS estimation is predominantly influenced by the assumed size of the observational uncertainties.

We emphasize that the observational uncertainties were introduced into the analysis by directly subjecting only the test dataset to an assumed level of “noise”, and the trained DNN model does not have prior knowledge of similar uncertainties in the training data. Nevertheless, these findings demonstrate that the DNN is able to accurately reconstruct the EOS from moderately noisy NS $M\left(R\right)$ observational data. The performance of the neural network can be further improved by introducing the measurement errors also in the training dataset. Since the systematic study of the measurement uncertainty effects is not the major focus of this work, detailed investigations, and also studying the effect of introducing “noise” in the training data, are left for following articles.

3.2. Application to Realistic EOSs

To test further the performance of the trained EOS DNN model, we apply it to several realistic EOSs from the CompOSE repository [154] (https://compose.obspm.fr (accessed on 25 March 2023)). CompOSE is an online tool that provides data tables containing state-of-the-art EOSs that can be readily used for various applications in astrophysics and nuclear physics. For the purpose of our analysis, we chose several EOSs within the range of the parameter space of the DNN training dataset: APR [134], BL [155], QMC-RMF2 [156], QMC-RMF3, [156], SK255 [157], and SK272 [157]. We generated the required input data following the procedure outlined in Section 2.3.1. For each EOS (Figure 7, left window), we calculated the $M-R$ relation (Figure 7, right window), and then drew ${10}^5$ random ensembles of 50 $(M_i,R_i)$ points to determine the mean and MAE of the reconstructed EOSs. The results of this test are shown in Figure 8 and demonstrate the ability of the EOS DNN to reconstruct realistic EOSs from $M\left(R\right)$ data. In all frames, the solid blue lines denote the ground-truth EOS and the red dot characters represent the mean of the DNN predictions respectively. The error bars represent the MAE with which the EOS is reconstructed in each case due to the random drawing of the $M\left(R\right)$ data, and in order to clearly separate the error contribution of this effect alone, they do not include the effect of assumed observational uncertainties. Among the realistic EOSs we considered, it is seen that the predicted $P\left(\rho \right)$ relation matches almost exactly the ground-truth values for the APR, BL, QMC-RMF2 and QMC-RMF3 models, while the the SK255 and SK272 EOSs are reconstructed less accurately, but still within the reconstruction errors.

Here, we briefly recall the main features of the realistic EOSs used in our analysis. The APR EOS was calculated using variational approaches with the A18 + delta v + UIX* interaction [134]. The BL equation of state (EOS) was derived from two-body and three-body nuclear interactions using $\chi$EFT, which includes the $\Delta \left(1232\right)$ isobar intermediate state [155]. The EOS has been derived using the Brueckner–Bethe–Goldstone quantum many-body theory in the Brueckner–Hartree–Fock (BHF) approximation with the continuous choice for the auxiliary single particle potential. The QMC-RMF2 and QMC-RMF3 equations of state (EOSs) were obtained using a relativistic mean-field (RMF) theory constrained by $\chi$EFT calculations of pure neutron matter within the density range of 0.08 fm${}^{-3}$ to 0.32 fm${}^{-3}$, as well as by the properties of isospin-symmetric nuclear matter around ${\rho }_0$ [156]. The SK255 and SK272 EOSs are unified models by Gulminelli and Raduta [157], computed with the SK255 and SK272 effective interactions [158]. The APR and BL EOSs are microscopic, while the rest of the EOSs are based on phenomenological models.

These results clearly demonstrate that a DNN, trained on a relatively simple dataset generated with the EOS metamodel discussed in Section 2.1, is able to generalize the task of reconstructing the EOS and predict accurately realistic EOSs.

3.3. Deducing Nuclear Matter Properties

3.3.1. Performance on the Test Dataset

In the following analysis, we examined the effectiveness of the trained NuPRO DNN model in extracting selected nuclear matter properties from the EOS of $\beta$-stable neutron-star matter. After the model was trained and the optimal architecture was determined (as shown in Table 2), we assessed its final performance by evaluating it on a test dataset composed of 5000 samples of $P\left(\rho \right)$, and corresponding sets of selected nuclear matter properties matching each EOS sample. To evaluate the model’s performance, we computed the standard deviation, ${\sigma }_{{ϵ}_i}$, of the residuals, ${ϵ}_i=Q_i^{DNN}-Q_i$, for each of the nuclear matter parameters: $K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$. Here, $Q_i$ is one of the selected five nuclear matter properties. By examining the standard deviation, we can determine the degree of accuracy and precision of the model’s predictions for each of these parameters.

In Figure 9, we present the results of our evaluation of the performance of the trained NuPRO DNN model in extracting selected nuclear matter properties from the EOS of $\beta$-equilibrium NS matter. We provide scatter plots of the distribution of the residuals for each of the nuclear matter parameters, along with the numerical values for the mean ${\mu }_{ϵ}$ and the standard deviation ${\sigma }_{ϵ}$. These results clearly demonstrate that the trained NuPRO DNN model achieved a high degree of accuracy in extracting the nuclear matter parameters. We noticed that the lower-order terms were extracted more accurately than the higher-order terms, and this is because the lower-order terms have a smaller range of possible values. The smaller range of possible values allows for a better interpolation precision of the model, which in turn leads to more accurate extraction of the lower-order terms. These results are also summarized in

Table 3. Mean, ${\mu }_{{ϵ}_i}$, and standard deviation, ${\sigma }_{{ϵ}_i}$, of the residuals, ${ϵ}_i=Q_i^{DNN}-Q_i$, of the trained NuPRO DNN (nuclear matter properties deep neural network) model, determined on the test dataset. All values are given in MeV. See text for details.

${\mathit{Q}}_{\mathit{i}}$	${\overline{\mathit{ϵ}}}_{\mathit{i}}$	${\mathit{\sigma }}_{{\mathit{ϵ}}_{\mathit{i}}}$
$K_0$	−0.22	2.77
$J_0$	0.79	4.09
L	0.23	0.49
$K_{sym}$	0.67	3.37
$J_{sym}$	0.57	5.34

.

3.3.2. Reconstructing $E_{sym}\left(\rho \right)$

Next, we focus on evaluating the ability of the trained NuPRO DNN model to accurately extract the nuclear symmetry energy parameters, namely L, $K_{sym}$, and $J_{sym}$, which are used to reconstruct $E_{sym}\left(\rho \right)$. The dense matter EOS is highly dependent upon the nuclear symmetry energy, which is a vital but uncertain component, and therefore it is imperative to explore whether DL techniques could offer a viable means of deducing it from astrophysical observations of neutron stars.

In our previous work [121], we pioneered the use of DL methods for extracting the nuclear symmetry energy from a set of neutron star observations in the $M-R$ or $M-\Lambda$ planes. In our “proof-of-concept” study [121], our main focus was on the extraction of $E_{sym}\left(\rho \right)$, and thus we generated our datasets by holding all parameters in Equation (2), representing the energy of symmetric nuclear matter, constant, and varying only the nuclear symmetry energy parameters L and $K_{sym}$ in Equation (3). We demonstrated that, under the given model assumptions, DNNs could extract $E_{sym}\left(\rho \right)$ effectively and accurately directly from astronomical observations of neutron stars. In the present investigation, we have advanced our deep learning methodology for extracting $E_{sym}\left(\rho \right)$ by significantly enlarging the parameter space of our neural network training dataset. To achieve this, we have kept only the parameters $E_0$ and $S_0$ fixed at their maximum likelihood values while varying the other parameters, namely $K_0$, $J_0$, L, $K_{sym}$, and $J_{sym}$, in Equations (2) and (3), to generate numerous samples of $P\left(\rho \right)$. As depicted in Figure 7, the augmented parameter space of the neural network training datasets also permits the modeling of predictions of modern realistic equations of state, which satisfy the constraints from recent mass–radius observations of neutron stars.

Figure 10 displays results for five instances selected from the test dataset. The reconstructed nuclear symmetry energy (represented by broken colored lines) for each input $P\left(\rho \right)$ sample closely matches the true $E_{sym}\left(\rho \right)$ (indicated by solid black lines). Similar outcomes were obtained for the remaining test data samples. We emphasize that the reconstructed nuclear symmetry energy is deduced by substituting the estimated values of the nuclear symmetry energy parameters, namely L, $K_{sym}$, and $J_{sym}$, predicted by the NuPRO DNN, into Equation (3). Moreover, we assumed that the $\beta$-equilibrium NS EOS, $P\left(\rho \right)$, is already known with a certain level of precision, such as being extracted from mass–radius observations of neutron stars by using the EOS DNN trained model.

3.3.3. Model Uncertainty

Let us not forget that realistic observations of neutron stars unavoidably harbor uncertainties, which in turn give rise to uncertainties in the inferred EOS of $\beta$-stable NS matter, and consequently in the extracted nuclear matter parameters and the symmetry energy $E_{sym}\left(\rho \right)$. To investigate the impact of errors in the reconstruction of the EOS on the inferred nuclear matter parameters and symmetry energy, we incorporated “noise” into the $P\left(\rho \right)$ data samples that portray the EOS and evaluated the nuclear matter parameters and $E_{sym}\left(\rho \right)$. We conducted experiments by varying the level of noise and investigated the resulting effect on the accuracy of the extracted nuclear matter parameters and symmetry energy. As shown in Figure 11, we elucidate the effect of introducing a 20% uncertainty to the input $P\left(\rho \right)$ data samples on the reconstructed $E_{sym}\left(\rho \right)$. In the left window, we show the exact EOS (represented by the solid blue line) and an EOS data sample containing 20% uncertainty (indicated by the red broken line). The reddish colored band denotes the uncertainty of $P\left(\rho \right)$. In order to assess the uncertainty in determining the nuclear matter parameters and the symmetry energy, we generate ${10}^5$ random sets of 50 equally spaced points in $\rho$, $P\left({\rho }_i\right)$, where $i=1,2,\dots ,50$, lying within the uncertainty band. We subsequently compute the mean and standard deviation for each of the nuclear matter parameters for each set. Thereafter, with every estimated set of nuclear matter parameters, we determined $E_{sym}\left(\rho \right)$ through Equation (3). The reconstructed symmetry energy is illustrated in the right window of Figure 11. The reddish colored band represents the mean absolute error (MAE) in deducing $E_{sym}\left(\rho \right)$, the solid line depicts the exact symmetry energy, and the red broken line indicates the mean symmetry energy. As anticipated, since the inferred symmetry energy is reconstructed via Equation (3), it closely follows the exact $E_{sym}\left(\rho \right)$ in a qualitative manner. Quantitatively, the estimated values begin to deviate moderately from the exact ones at approximately $\rho \ge 2{\rho }_0$; however, they remain within the range specified by the mean absolute errors of the model for the assumed uncertainty of the input EOS. The mean, standard deviation, and MAE for each of the nuclear matter parameters for the specific example shown in Figure 11 are presented in

Table 4. Values for the exact, predicted, mean ${\mu }_i$, standard deviation ${\sigma }_i$, and mean absolute error (MAE) of the nuclear matter parameters $Q_i$ = [$K_0,J_0,L,K_{sym},J_{sym}$] for the example illustrated in Figure 11. All values are given in MeV. Please see the text for further details.

${\mathit{Q}}_{\mathit{i}}$	Exact	Predicted	${\mathit{\mu }}_{\mathit{i}}$	${\mathit{\sigma }}_{\mathit{i}}$	MAE
$K_0$	259.19	256.49	292.18	55.66	50.50
$J_0$	−78.38	−73.22	−69.91	69.05	56.66
L	58.90	58.99	71.38	12.79	13.35
$K_{sym}$	−225.75	−223.30	−227.52	51.41	42.23
$J_{sym}$	520.61	516.33	531.51	226.82	179.65

. The results are highly analogous for the remainder of the data samples from our test dataset.

It is important to note that the uncertainties presented in our analysis were solely introduced to the test dataset, and the trained DNN model does not possess any prior knowledge of uncertainties in the training data. Despite this, our findings demonstrate that the neural network is capable of accurately extracting the nuclear matter parameters and reconstructing $E_{sym}\left(\rho \right)$, even when faced with moderately noisy input EOS data. In order to further improve the performance of the DNN, it may be beneficial to introduce uncertainties to the training dataset as well. Although the impact of measurement uncertainties was not the primary focus of the current study, the potential effects of such uncertainties can be studied in detail in future works. This includes investigating the effects of introducing “noise” to the training data, as well as conducting systematic studies on the impacts of measurement uncertainties. These investigations may provide further insight into the behavior of the DNN and help to enhance its performance in future applications.

3.3.4. Application to Realistic Nuclear Models

So far we have demonstrated that the trained NuPRO DNN model performs with high accuracy on the test dataset. Having established this, we now proceed to applying the model to a set of realistic EOSs, which were previously discussed in Section 3.2. It is important to acknowledge the complexity of the inference task and limitations of our model assumptions before discussing the results. Firstly, the DNN model was trained on a dataset that assumes the nuclear matter EOS depends on the matter isospin asymmetry via a quadratic dependence only, as described in Equation (1). Secondly, in regards to the hadronic part of the EOS, we employed the models parameterized by Equations (2) and (3) for symmetric nuclear matter and nuclear symmetry energy $E_{sym}\left(\rho \right)$, respectively. These models were used in the density range of approximately 0.04 fm${}^{-3}$ to 0.8 fm${}^{-3}$. It is important to note that beyond the saturation density ${\rho }_0$, these expressions should be regarded solely as parameterizations, and not as Taylor expansions. Thirdly, our model was based on the assumption that $\beta$-equilibrated neutron-star matter consists of nucleons, electrons, and muons. This assumption was made to simplify the problem, and it may not accurately represent the composition of matter in neutron stars, which may include other exotic particles, such as hyperons, or quark matter.

In Figure 12, we present the discrepancies between the predicted values of the nuclear matter parameters $K_0$, L, and $K_{sym}$ obtained from the trained NuPRO DNN model and their true values obtained from the CompOSE repository for each of the considered realistic EOSs. These discrepancies represent the residuals in our analysis. We note that $K_{sym}$ values were not available for the APR and BL EOSs, and hence, we do not show the residuals for these models. We also computed the standard deviations, ${\sigma }_i$, of the residuals to evaluate the uncertainty in estimating the nuclear matter parameters from a realistic neutron-star EOS in $\beta$-equilibrium. The standard deviations for $K_0$, L, and $K_{sym}$ are 30.18 MeV, 11.22 MeV, and 19.09 MeV, respectively. These values are smaller than the reported uncertainties in the literature [159].

Here we reiterate that our model assumptions and limitations should be taken into account when interpreting these results. Therefore, the applicability of our results to other types of matter, such as hyperonic matter or quark matter, remains an open question. Nonetheless, our findings demonstrate the potential of using DNN models to extract nuclear matter parameters from astrophysical observations of neutron stars.

Precise measurements of the masses and radii of a sufficient number of neutron stars would ultimately allow for the accurate determination of the EOS of $\beta$-stable matter through converting the $M\left(R\right)$ curve, via various methods, to the underlying EOS [119]. However, extracting the properties of nuclear matter from the $\beta$-equilibrated equation of state presents another challenge since the exact composition of a neutron star’s interior is unknown, and even determining the proton fraction is highly challenging [119]. For instance, in Ref. [160] the authors applied a Bayesian approach and showed that it is impossible to determine precisely the properties of nuclear matter from the EOS of $\beta$-stable matter. Similarly, the authors of Ref. [161] demonstrated the existence of multiple solutions for the determination of the NS interior composition from the $\beta$-stable matter EOS, owing to the high level of degeneracy. Additionally, accurately determining the nuclear symmetry energy from the $\beta$-equilibrium EOS presents a challenge as it necessitates precise knowledge of the equation of state of symmetric nuclear matter [162]. This is because determining the proton fraction in the interior of a neutron star is dependent on an accurate understanding of symmetric nuclear matter. Due to the inherent correlations among the nuclear matter properties and their intertwined relationship to both the EOS of symmetric nuclear matter and the nuclear symmetry energy, the precise determination of the nuclear matter properties from the beta-stable matter EOS is extremely challenging [163]. In order to determine unambiguously the nuclear matter properties from the EOS of beta-stable matter, both the EOS of symmetric nuclear matter and the nuclear symmetry energy need to be known independently at a given density [163]. These considerations underscore the key role and potential of the data-driven methods presented in this work as they provide a model-independent avenue to deducing the EOS of $\beta$-stable matter, and in turn, the nuclear matter properties and $E_{sym}\left(\rho \right)$.

4. Summary and Outlook

In this work, we have presented a successful application of deep learning techniques to extract the equation of state of dense matter directly from neutron star observations. By analyzing simulated mass and radius measurements of neutron stars, our results demonstrate that deep neural networks can accurately extract the EOS of $\beta$-stable matter. Moreover, the findings presented in this work highlight the capability of the trained DNN model to determine selected nuclear matter properties, including the nuclear symmetry energy $E_{sym}\left(\rho \right)$. Most importantly, we have demonstrated that our DL approach is highly effective in accurately extracting realistic EOSs and associated nuclear matter properties from observational data of neutron stars. This achievement is significant, as it represents a major step towards the ultimate objective of determining the EOS of dense nuclear matter. Our results underscore the potential of DL-based techniques in the era of multi-messenger astrophysics, where a growing body of NS observational data is rapidly emerging.

In the near future, we plan to systematically examine the uncertainties associated with the NS observational data and the DNN model, and their impact on the model’s performance. By understanding the effect of these uncertainties, we aim to explore potential approaches to further enhance the model’s reliability and performance. In particular, in order to apply our approach practically, it is essential to consider the empirical errors and uncertainties and incorporate them consistently into the formalism. A possible strategy to achieve this is to recast the regression problem of extracting the EOS and nuclear matter properties into a statistical model.

In particular, in future studies, our aim is to use Bayesian neural networks to carry out the inference procedure. In this approach, rather than obtaining deterministic values, the weights of the network are described by probability distributions by applying a prior over the network weights [164]. In following works, we also plan to apply our DL approach to real observational data of neutron stars, which would enable the extraction of a model-independent EOS, nuclear matter properties, and symmetry energy. Lastly, we aim to explore likelihood-free inference approaches through the use of normalizing flows [165]. These methods have the ability to model intricate posteriors by applying non-linear transformations to a simple posterior shape, such as a multivariate Gaussian, without the need for direct likelihood evaluation.This approach has already generated considerable interest in the scientific community and it has been successfully applied in multiple research domains. For example, a recent research work [166] employed a likelihood-free inference technique that used normalizing flows to efficiently estimate the parameters of eight binary black hole (BBH) events from the first LIGO Gravitational Wave Transient Catalog, GWTC-1 [167]. As the next-generation space telescopes and GW detectors become operational, they are expected to detect and observe a large number of compact binary collisions and neutron stars spanning the history of the universe. These instruments will identify more than a million events each year, including thousands of BNS and NSBH detections, and processing this vast amount of observational data with speed and accuracy will be critical. It should be noted that traditional Bayesian inference methods may not be feasible for analyzing such large datasets and, hence, more modern techniques such as normalizing flow models may be essential in accurately and rapidly extracting important neutron star parameters.

Ultimately, as the number of observed neutron star events increases, these modern data-driven techniques will allow us to efficiently handle the growing volume of observational data and precisely determine the equation of state of dense nuclear matter and the nuclear symmetry energy.

Software: Matplotlib (https://matplotlib.org (accessed on 25 March 2023)), Seaborn (https://seaborn.pydata.org (accessed on 25 March 2023)), NumPy (https://numpy.org (accessed on 25 March 2023)), SciPy (https://scipy.org (accessed on 25 March 2023)), TensorFlow (https://www.tensorflow.org (accessed on 25 March 2023)), Adobe Illustrator (https://www.adobe.com/products/illustrator.html (accessed on 25 March 2023)).