SIGS - Seismic Inferences for Glitches in Stars

The increased amount of high precision seismic data for solar-like stars calls for the existence of tools that can extract information from such data. In the case of the study of acoustic glitches there are no publicly available tools and most existing ones require a deep knowledge of their implementation. In this work a tool is presented that aims to both simplify the interaction with the user and also be capable of working automatically to determine properties of acoustic glitches from seismic data of solar-like stars. This tool is shown to work with both the Sun and other solar analogs but also shows that are still severe limitations to the methods used, when considering smaller datasets.


Introduction
Acoustic glitches, such as the base of the convective zone (therein BCZ) and the helium second ionization zone (therein HeII), are regions in the interior of a star where the sound speed suffers an abrupt variation due to the sharp change in the internal structure. Since the eigenfrequencies of a star are directly related to the local sound speed, these glitches cause perturbations in them. These perturbations exhibit an oscillatory behavior when considering frequencies of sequential radial order [1][2][3]. It is thus possible to characterize some properties of these acoustic glitches by isolating and studying their contribution to the eigenfrequencies. This oscillatory signature from the acoustic glitches has been extensively studied for the Sun to determine the location of the two previosly mentioned acoustic glitches, BCZ and HeII [4,5] and it has also been applied to solarlike stars [6][7][8][9].
However, the tools existent nowadays to perform this study are not publicly available and require an intimate knowledge of the implemented methods in order to achieve good results. The work presented here aims to improve upon existing methods in order to make available a numerical procedure that is both more precise than the current methods and at the same time highly user-friendly so that it can be openly distributed to the community. Starting from the work developed by Faria [8], which in turn was adapted from the works of Monteiro et al. [4], Monteiro and Thompson [5] and Mazumdar et al. [9] the method was improved upon and an automatic pipeline was implemented to achieve robust results independent of user intervention.
The presented tool is validated initially using observational data from the Sun and afterwards it is applied to both stars from the 16 Cygni binary and also to 10 solar-like stars from chosen from the sample of Mazumdar et al. [9].
In the next section both methods implemented in the tool are described. Then in section 3 the code is validated using data from the Sun. In section 4 the automatic implementations of the code are described. In section 5 the results obtained for the mentioned stars are shown and finally in section 6 the results are discussed and some conclusions are drawn.

Methods
The program presented in this work is comprised of two similar but independent methods capable of isolating the oscillatory signature from the acoustic glitches and then determining some of their properties.

Method A
Method A isolates the signature in the frequencies themselves by iteratively removing a smooth component from the frequencies and then fitting the residuals to the free parameters from equation +A HeII ν r ν sin 2 (2πβ HeII ν)cos(4πτ HeII ν + 2φ HeII ).
(1) adapted from Faria [8] where A BCZ , τ BCZ and φ BCZ correspond to the amplitude, acoustic depth and phase of the signal from the base of the convective zone and A HeII , τ HeII , φ HeII and β HeII correspond again to amplitude, acoustic depth and phase and finally, the acoustic width of the signal, this time due to the helium second ionization zone. ν s is the smooth component removed from the frequencies which is obtained independently for the N l frequencies of each angular degree l by fitting a polynomial of degree N l − 1 with third derivative smoothing [4]. This smoothing can be controlled by defining a smoothing parameter λ which can be controlled by the user but is by default determined automatically by the code. After removing the smooth function from the frequencies, the residuals are fitted to Equation 1 using the PIKAIA genetic algorithm [10] which is a global minimization algorithm that utilizes concepts inspired by the processes of evolution by natural selection. This method is advantageous compared to hill climbing methods because it doesn't require an initial guess and instead it explores the user-defined range of possible values for each of the free parameters and should, given enough time, converge to global minimum in the parameter space. Minimization is achieved by initiating a population of sets of possible solutions chosen randomly from the interval and then evolving then for a specified number of generations. Through the generations the best sets from the population are chosen to move forward by using a standard chi squared fitness function and after the evolution process the best set of parameters is chosen as the correct one.
The process of executing the smooth removal and fitting the residuals to the equation is repeated until the smooth function is perfectly characterized and convergence is said to be achieved. Figure 1 shows a diagram of the general workflow of this method Figure 1. Diagram of the workflow of method A of SIGS which determines the acoustic depths of the base of the convective zone and helium II ionization zone of a star by extracting the signal caused by these glitches directly from the oscillation frequencies of the star.

Method B
In method B, instead of using the frequencies, the method searches for the signature of the glitches in the second differences [3].
In this method the smooth component is only removed once and the functional form of this component is a polynomial of a degree up to three depending on the amount of data available for the fit.
After removing the smooth function the procedure is similar to the other method with the exception that there is no iteration on the smooth component and the functional form of the signal is taken to be which is very similar to the one from the frequencies with the same seven free parameters but a different smooth function as mentioned previously. This expression was adapted from the work of Mazumdar et al. [9], which in turn was adapted from the work of Houdek and Gough [7] Figure 2 shows the workflow of method B.

Monte Carlo Simulations
To add robustness to both methods and estimate the errors in the measurements of the parameters we perform Monte Carlo simulations. We generate alternative data sets by taking the original frequencies and changing them by a value sampled from a normal distribution with standard deviation equal to the uncertainty of the original frequency.
The methods are then used on all the new data sets to obtain a distribution of values for the final parameters.

Validation
To validate the implementations of both methods we apply them to both low-degree solar frequencies [11] and then compare the results with that of previous works.
The results obtained for both methods are summarized in Table 1. Figures 3 and 4 show fits to the best parameters   The results obtained are successful in showing that the code is working correctly. The small difference in values compared to the literature is to be expected considering the small changes made to the methods implemented and to the more recent data used.

Automatic Pipeline
After validating the code using the Sun, the program was used to determine the parameters of a grid of models and the results were used to build a pipeline that was capable of reducing the parameter interval automatically to improve robustness to the results without requiring user interaction.
A grid [12] of CESAM [13] models was considered with masses between 0.8M and 1.2M in intervals of 0.04M and each mass was considered with 5 different ages (1.5, 3, 4.5, 6 and 8 Gyr). The frequencies of the models were computed using the posc algorithm [14]. The results from this run were then used to associate the obtained values of the acoustic depths of both glitches to the values of effective temperature (Teff) and large frequency separation (Large_sep) of the stars through a multivariate quadratic regression. The result of the regression is shown in Figure 5 . Only the results using method A were shown since the results for method B are very similar. The results from this regression allowed for the definition of a function that estimates the value of the acoustic depths of both glitches from the Teff and Large_Sep of a star, improving the focus of the fitting procedure.

Results
To test the complete tool which is to be available to the community the program was applied to the stars from the 16 Cygni binary, with data from Verma et al. [15] and to 10 solar-like stars from the sample of Mazumdar et al. [9] with frequencies from Appourchaux et al. [16]. The errors determined for the results were obtained by executing the program for 500 realizations of the observational data of each star, using the aforementioned Monte Carlo simulations. The results for 16 Cyg A and B are presented in Table 2. Table 3 shows the results for the 10 selected stars.

Discussion and Conclusions
Regarding the stars from the 16 Cygni binary, both show precise results that are very similar with the ones found in other works [15]. However, for the 10 selected stars the results are considerably poorer with some stars having no convergence achieved with any of the methods. Nevertheless, the results are very similar to those found in Mazumdar et al. [9] which might indicate that the poor results are related to the amount of available data for these 10 stars and not to the methods adopted, since these are proven to work for stars with more data. To conclude, from the results obtained it seems that, even with the automatic pipeline that reduces the parameter space for the fitting algorithm, the methods are highly dependent on the amount and quality of the available data. Still, the results for the Sun and the stars of the 6 Cygni binary show that is tool has potential as it achieves the same results as other methods without the need for user intervention and guidance. There are also some improvements that can be done to the program such as increasing the number of free parameters to the fitted functions to try and constrain the signal better.