We trained denoiser autoencoding neural networks on medium resolution simulated optical spectra of latetype stars to demonstrate that the reconstruction of the original flux is possible at a typical relative error of a fraction of a percent down to a typical signaltonoise ratio of 10 per pixel. We show that relatively simple networks are capable of learning the characteristics of stellar spectra while still flexible enough to adapt to different values of extinction and fluxing imperfections that modifies the overall shape of the continuum as well as to different values of Doppler shift. Denoised spectra can be used to find initial values for traditional stellar template fitting algorithms and  since evaluation of pretrained neural networks is significantly faster than traditional template fitting  denoiser networks can be useful when a fast analysis of the noisy spectrum is necessary for example during observations between individual exposures.
['Astrophysics - Instrumentation and Methods for Astrophysics', 'Astrophysics - Solar and Stellar Astrophysics']
Denoising medium resolution stellar spectra with neural networks
{'No Header': 'DOI: xxx/xxxx', 'Denoising medium resolution stellar spectra with neural networks': 'Balázs Pál* 1,2,3 | László Dobos 3  \n- 1 Department of Physics of Complex Systems, ELTE Eötvös Loránd University, Hungary\n- 2 Institute for Particle and Nuclear Physics, HUN-REN Wigner Research Centre for Physics, Hungary  \n3 Department of Physics & Astronomy, The Johns Hopkins University, Maryland, USA', 'Correspondence': '*Balázs Pál, Department of Physics of Complex Systems, Eötvös Loránd University, Pázmány Péter avenue 1/A, Budapest, Hungary. Email: pal.balazs@ttk.elte.hu  \nWe trained denoiser autoencoding neural networks on medium resolution simulated optical spectra of late-type stars to demonstrate that the reconstruction of the original flux is possible at a typical relative error of a fraction of a percent down to a typical signal-to-noise ratio of 10 per pixel. We show that relatively simple networks are capable of learning the characteristics of stellar spectra while still flexible enough to adapt to different values of extinction and fluxing imperfections that modifies the overall shape of the continuum, as well as to different values of Doppler shift. Denoised spectra can be used to find initial values for traditional stellar template fitting algorithms and - since evaluation of pre-trained neural networks is significantly faster than traditional template fitting - denoiser networks can be useful when a fast analysis of the noisy spectrum is necessary, for example during observations, between individual exposures.', 'KEYWORDS:': 'stars: atmospheres, methods: numerical, methods: data analysis, techniques: spectroscopic', '1 INTRODUCTION': 'Sky surveys focusing on observations of resolved stars often have to deal with the apparent brightness of target objects covering an interval of 6-7 magnitudes. Fiber-fed spectrographs on 8-meter-class telescopes can observe thousands of targets simultaneously, with the same exposure time, resulting in a broad range of signal-to-noise ratios for the final spectra. Often, the low surface density of bright stellar targets, combined with the limitations of fiber placement and the available telescope time, balances the target lists toward fainter magnitudes. This poses a challenge in both flux calibration and physical parameter inference. Atmospheric parameters, as well as the line-of-sight velocity of stars are determined by fitting theoretical stellar spectrum templates to the observed spectra which fitting process is prone to failure due to low signal-to-noise or incorrect flux calibration.  \nThe main motivation to our work is the demanding requirements of Galactic Archaeology (GA) which investigates the structure and formation as well as the dynamical and chemical evolution of the Milky Way, the Andromeda Galaxy and their  \nsurrounding satellite galaxies. Low and medium resolution spectroscopy is an important tool to gather this information about resolved stars and new generation instruments such as DESI (DESI Collaboration et al., 2022; Dey et al., 2023) and Subaru PFS (Sugai et al., 2014; Takada et al., 2014) are capable of observing red-giant branch (RGB) stars as far as the Andromeda Galaxy (Dey et al., 2023) and even main sequence turn-off stars at the distance of some Milky Way satellites. By gaining insight into the dynamics of these systems via measuring the line-of-sight velocities of individual member stars, and understanding the chemical evolution by analyzing elemental abundances, we can learn about the small scale characteristics of the dark matter and its early time interaction with ordinary matter (Takada et al., 2014).  \nThe low brightness of the targets - with the sample median easily as faint as 22 . 5 magnitudes in the 𝑖 -band, which translates to a typical 𝑆 / 𝑁 of 5 per pixel or 10 per resolution element at medium resolution with 3 hours of exposure time - poses a significant challenge for the analysis of these observations. Most stellar spectrum analysis methods rely on a continuum normalization step which first fits a template or analytic continuum to the observed flux and then divides the flux  \nwith the fitted continuum to isolate the spectral lines. Classical methods often struggle to find the continuum and best fit stellar template to noisy observations which can lead to non-converging optimizations and strongly biased parameter estimates.  \nIn this paper, we investigate if deep neural networks, more precisely denoising autoencoders (DAEs), are capable of recovering the underlying theoretical synthetic models from simulated medium resolution optical spectra with different levels of noise. We train and evaluate our DAE on a training set consisting of a 100 , 000 simulated stellar spectra of varying atmospheric parameters and 𝑆 / 𝑁 . Once neural networks are trained, they can be evaluated on an observation in a few milliseconds, hence denoising autoencoders can act as a fast, on-the-fly algorithm to find the approximate best-fit synthetic spectrum and supplement more traditional algorithms to find good initial values for iterative optimization algorithms, hence reduce computation time.', '1.1 Denoising autoencoders': "Autoencoders (Kramer, 1991) are neural networks that have the same, usually large number of inputs and outputs that are connected with hidden layers that have an hourglass shape with a hidden middle layer of only a small number of nodes, referred to as the bottleneck layer . During training an autoencoder, the requirements is that it reproduces its input on the output layer as close as possible (Hinton & Salakhutdinov, 2006). The first half of the autoencoder, between the inputs and the middle layer is referred to as the encoder and has a narrowing shape where each layer is followed by a layer with a smaller number of neurons. This first half of the network is responsible for encoding the input into as few numbers as the number of neurons in the bottleneck layer. The encoder is a non-linear mapping between the high-dimensional input and a low dimensional latent space . The output of the bottleneck layer is a lossy compression of the input, where the compression depends on the information content of the input, the network architecture and how well the network is trained. The second half of the network, the decoder , learns how to reconstruct the original input from the compressed representation of the bottleneck layer and it usually consists of layers with an increasing number neurons toward the output.  \nDenoising autoencoders are very similar in concept to regular autoencoders (AEs) but are trained differently and can show differences in architecture (Vincent, Larochelle, Bengio, & Manzagol, 2008). In case of a denoising network, the low information flux through the bottleneck layer is primarily responsible for the denoising capabilities.  \nThe weights of the autoencoders are trained by requiring that the network copies the input data to its output as accurately  \nTABLE 1 Summary of the observational parameters and parameter ranges that were sampled uniformly to generate the simulated spectra. These parameters are very similar to the typical conditions at at Mauna Kea Observatory, except for the seeing which tends to be 0 . 8 '' or higher. The exposure time was set to 5 hours, with 20 exposures and detector read-outs, 15 minutes each.  \nas possible. In contrast, denoising autoencoders are trained by requiring that the network produces a noiseless, or at least less noisy output from a noisy input. Since the primary goal of DAEs is to eliminate noise as much as possible, rather than encoding the input into a compressed representation as traditional AEs, the network architecture of DAEs does not necessary follow the hourglass shape and the bottleneck layers usually contains more neurons as in case of AEs.", '2 SIMULATING STELLAR SPECTROSCOPIC OBSERVATIONS': 'To train the denoising autoencoder, we generated a large ensemble of simulated spectra based on the BOSZ synthetic spectrum grid (Bohlin et al., 2017). The observational parameters were chosen to match the typical environmental conditions at the Subaru Telescope at the Mauna Kea Observatory on dark and gray nights, as summarized in Tab. 1 . Each of the parameters is sampled uniformly within the specified intervals.  \nWe simulated the spectra as if they were observed with an instrument very similar to the 𝑅 ≈ 5000 medium resolution mode of the Subaru Prime Focus Spectrograph (Tamura et al., 2016). The parameters of the spectrograph are listed in Tab. 2 . The instrument was specifically designed for Galactic Archeology and its wavelength range was chosen to cover lines that carry the most information about line-of-sight velocity and elemental abundances in a limited wavelength range.  \nTABLE 2 Summary of the spectrograph parameters similar to the medium resolution of the Prime Focus Spectrograph of the Subaru Telescope.  \nTABLE 3 Summary of the of the fundamental stellar atmospheric parameters and other simulation input parameters. The intervals are sampled uniformly in their respective ranges to generate the simulated spectra. The fundamental parameters cover the range of interest of Galactic Archeological observations but also include non-physical models that we did not exclude.  \nThe ranges of fundamental stellar parameters were chosen to cover stellar types from M to G, which differ from the typical targets of Galactic Archeology, which prefers the stellar types Gand F. We chose the lower temperature range, because of the higher variance of the absorption features, which makes the denoising task more challenging and exercises the capabilities of the neural network to learn the spectral features.  \nSince the parameters are uniformly and independently drawn from continuous intervals, the simulated spectra include many non-physical models which should have no significant impact on the performance of the neural networks trained with them. The input parameters of the simulations are summarized in Tab. 3 . We used these values to generate both, the training set and the validation set for optimizing the weights of the denoising autoencoder, as detailed in Sec. 3.2.  \nThe spectrum simulation process can be separated into two main parts. The first part is the interpolation of the synthetic stellar grid to calculate a model spectrum at a randomly chosen  \nset of atmospheric parameters and the second part is the simulation of all observational and instrumental phenomena that affect the light of the star between leaving its atmosphere and getting detected by the instrument. The output of the simulation is a noiseless spectrum with the specific flux expressed in physical units, as well the variance of the flux for each spectral pixel. The variance vector is later used to generate realizations of the noise during the training process of the neural network.  \nIn addition to the precomputed observational effect, we apply further transformations to the simulated spectra during training of the neural network which serve the purpose of dataset augmentation. For instance, we generate new realizations of the photon noise on-the-fly, at each training epoch as well as introduce a wavelength-dependent systematic flux calibration error to mimic the effects of improper fluxing. These will be detailed in Sec. 3.3.  \nThe observational and instrumental effects are modelled according to the Exposure Time Calculator (ETC) of the Subaru PFS. The Subaru PFS ETC is an adaptation of the original WFIRSTETCcode, as detailed in Hirata et al. (2012). In order to be able to generate hundreds of thousand of spectra, we modified the ETC code to output separately (i) the atmospheric and instrumental transmission functions, (ii) the individual noise terms originating from the night sky, as well as (iii) the instrumental line spread function (LSF) as a function of wavelength. In turn, the LSF, at each wavelength pixel, was fitted with a Gaussian function with a wavelength-dependent width so that it can be evaluated at any wavelength and spectral resolution, consistently with the wavelength grid of the input synthetic spectra.  \nWe ran the modified ETC on a grid of observational parameters and ingested the outputs into data arrays. Our parallel Python implementation of the observation simulator interpolates these precomputed grids to calculate and sum up the noise terms that contribute to the flux variance.  \nBelow, we discuss the simulation steps in order, and explain further processing that happens during the training process only.', '2.1 Spectrum interpolation': 'We start the simulations by interpolating a regular grid of fluxed synthetic stellar spectra which are higher resolution than the instrument in order to be able to confidently evaluate the sensitivity of the denoiser to sub-pixel Doppler shifts. We chose to use the BOSZ synthetic stellar grid (Bohlin et al., 2017) of ATLAS9 models (Kurucz, 1979; Mészáros et al., 2012). The BOSZ grid is available at a broad range of resolutions and its underlying model and line lists are considered reliable in the parameter ranges we are interested in. We interpolate the flux using one-dimensional cubic splines in the  \ndirection of a single parameter only, chosen randomly from the four fundamental parameters [Fe/ 𝐻 ] , 𝑇 eff , log 𝑔 and [ 𝛼 /Fe] , meanwhile the rest of the parameters are kept confined to grid values.  \nWhen we later evaluate the performance of the denoising autoencoder, as a function of the each of fundamental parameters, we will only consider the models that are interpolated in the direction of the particular parameter.', '2.2 Observational and instrumental effects': 'After interpolating the synthetic stellar grid, we apply each observational effect in the specific order as follows.  \n- 1. Draw a value for the Doppler shift from a uniform distribution between -0 . 003 ≤ 𝑧 ≤ 0 . 003 , and shift the wavelength of the interpolated model accordingly, without rebinning or resampling the spectrum.\n- 2. Draw a value for the extinction randomly from a uniform distribution between 0 ≤ 𝐸 ( 𝐵 -𝑉 ) ≤ 0 . 1 and apply the mean extinction curve of the diffuse dust in the Milky Way, as determined by Cardelli, Clayton, and Mathis (1989), to the Doppler-shifted spectrum.\n- 3. Convolve the spectrum, still at its original high resolution, with a line spread function similar to the LSF of full optical system of Subaru PSF, including the pixelization of the detector.\n- 4. Interpolate the spectrum to the wavelengths of the detector pixels. In our case, the number of detector pixels is 4096 .\n- 5. Given a randomly drawn apparent magnitude, calculate the flux in physical units as well as the photon count from the source and the different components of the night sky including the continuum, sky lines and reflected sunlight from the Moon.\n- 6. Calculate the flux variance in each detector pixel from the photon counts and detector characteristics as if all exposures were stacked together into a single spectrum.\n- 7. Normalize the flux to have a median value of 0 . 5 in the wavelength range of 6500 -9500 Å .  \nFinally, we measure the signal-to-noise as the median of the ratio of the flux to the square root of its variance in each pixel over the entire wavelength coverage of the instrument. Instead of using 𝑆 / 𝑁 as the input parameter, we normalize the synthetic flux of the model spectrum to a magnitude drawn from a uniform distribution that covers the typical magnitude range of potential targets. Although this method will not result in uniform sampling in 𝑆 / 𝑁 , it will make sure that the flux error  \nFIGURE 1 An example of a noiseless simulated stellar spectrum observation and its noisy counterpart. A noise component was added to the original spectrum as described in Sec. 3.3. The spectrum is rescaled by setting the median of the flux in the 6500 -9500 Å wavelength range to 0 . 5 . The 𝑆 / 𝑁 value indicated in the bottom left corner is the median pixel-wise signal-to-noise ratio over all pixels.  \n<!-- image -->  \nfrom sky subtraction, scattered light and the noise contribution of the detector are taken into at account at the right ratio. Rescaling the noise level at a later step to a given 𝑆 / 𝑁 would change the ratio of the noise components.  \nReal observations contain many pixels that need to be masked due to various reasons such as cosmic rays, bright sky lines, low detector efficiency etc. In this paper, we do not mask pixels of the simulated spectra and do not train the denoising autoencoder to work with masked inputs.  \nThe final normalization step is necessary to transform the flux values of the simulated spectra into the [0; 1] interval, which our deep neural network with a sigmoid activation function at its output can handle.  \nWe assume that the multiple exposures of the same target were taken at the same observation conditions, hence we only run the simulation once and assume perfect stacking, with perfectly matching wavelengths and LSF when calculating the variance of the flux.  \nFig. 1 shows a typical simulated spectrum of a star with parameters [Fe/H] = -1 . 5 , 𝑇 eff = 4183 K, log 𝑔 = 2 . 0 , [ 𝛼 /Fe] = 0 . 0 , with and without a random realization of the noise added to the flux. The purpose of the denoising autoencoder is to recover the noiseless spectrum (blue curve) from the noisy output (gray curve).  \nIn addition to the observational effects detailed above, during the training process, as part of the data set augmentation, we multiply each spectrum with a different, randomly generated function that varies slowly with wavelength. This was done to simulate the effects of incorrect flux calibration. The function is generated using a random walk process, filtered with a low-pass filter, resulting in a higher-order function  \nthan typical polynomial corrections used for flux calibration. Although the flux calibration effects are most often systematic, training the network with a randomized flux calibration error will make it resilient against both random and systematic calibration errors. The maximum amplitude of the flux calibration error was set to be 2 per cent, which is a typical value (Bohlin, 2014).', '3 DENOISER IMPLEMENTATION': 'To realize a denoising autoencoder, we trained a relatively shallow, hourglass-shaped, fully connected network with rectified linear unit activation functions between the hidden layers, and a sigmoid function before the output using stochastic gradient descent with no momentum. The optimization of the weights was based on a training and a validation data set with 100 , 000 simulated spectra in each. In addition, we augmented the data sets with a different realization of the flux noise at each training step. During training, the gradient at each step was evaluated by setting the network input to the noisy spectrum, while the output was required to be the noiseless spectrum.  \nThe network architecture, data set augmentation and training process are described in detail in the following subsections.', '3.1 Neural network architecture': 'The architecture of the denoising autoencoder we used in this study is illustrated in Fig. 2 . Both, the encoder and the decoder side of the network consist of 4 -4 fully-connected layers, with an input and output size of 4096 , the number of wavelength pixels of the simulated spectroscopic observations. The number of neurons in each hidden layer of the encoder decreases linearly, while on the decoder side, the number of neurons increases the same way. Each fully-connected layer is followed by a bias vector with the same size as the output of that layer. After each layer, except prior to the output, a rectifier linear unit (ReLU) activation was used. Before the output, a smoother sigmoid activation was inserted instead, to produce the final output.  \nBefore training, the weights, as well as the bias vector of the fully-connected layers are initialized using the Xavier initialization scheme with values drawn from the uniform distribution  (-√ 𝑘, √ 𝑘 ) , where 𝑘 = 6/(inputs + outputs) , a function of the number of inputs and outputs of the layer. The size of the bottleneck layer was set to a relatively large number of 256 neurons.  \nIn summary, all layers in our autoencoder network follow the conventional fully connected feed-forward neural network  \narchitecture as described by  \n𝐙 = 𝜎 ( 𝐖𝐱 + 𝐛 ) , (1)  \nwhere 𝐱 is the input vector, 𝐖 is the weight matrix, 𝐛 is the bias vector, 𝜎 is the activation function, with 𝐙 being the output vector of the layer.  \nTo implement the network, we used the PyTorch deep learning framework (Paszke et al., 2019).', '3.2 Training and validation set properties': 'Using the stellar parameter distributions described in Sec. 2 and summarized in Tab. 3 , we generated separate data sets for training and validation, each with 100 , 000 simulated spectra. As we have used entirely synthetic data sets, we did not use a third, dedicated test set for the evaluation of the trained network, instead we used the validation set for this purpose.  \nAs mentioned in Sec. 2.2, the spectra, and the corresponding flux variance array were normalized, so that flux values cover the [0 , 1] interval, thus no further input or output normalization happens during the training of the network. However, we point out that the denoiser is unable to produce an absolutely noiseless spectrum. When some remaining noise is present that extends beyond the [0 , 1] interval (i.e., part of the predicted spectrum is outside this interval), the output values are clipped to this range. This is because the output of the network is bounded by the sigmoid function at the end of the network, which produces values in the [0 , 1] interval. Consequently, this results in a slightly non-Gaussian residual noise in the output.', '3.3 Data augmentation': "In addition to generating a large number of simulated spectra, we further augment the number of training samples by onthe-fly transformations of the input and the output during the training process. The steps of the augmentation process are as of follows.  \n- 1. Generate a random realization of the flux calibration error and apply it to both the input and output.\n- 2. Generate a random realization of the observational noise and add it to the input spectrum. The output is kept noiseless.\n- 3. Generate random values for the scale and offset and apply to both the input and the output.  \nAt each training epoch, we generate a new noise realization for every spectrum based on the flux variance vector. Using new realizations of the noise at each training epoch is a very efficient way of preventing overfitting the training set and,  \nFIGURE 2 Schematic overview of our denoising autoencoder and the training pipeline including the augmentation steps. The noiseless input spectra, labeled as 'original spectrum', were simulated as described in Sec. 2. During the training process, these spectra are fed into an augmenter, described in Sec. 3.3, resulting in the 'noisy measurement'. Following this, the noisy observations are processed by the DAE network, which is trained to output the denoised spectra, labeled here as 'denoised spectrum'. During the training, a stochastic gradient descent optimizer is used to find the network weights that minimize the RMSE loss calculated from the 'denoised spectrum' and the 'original spectrum'.  \n<!-- image -->  \nsince noise in inherent to spectroscopic observations, it is also a very natural method of data augmentation.  \nIn addition to adding noise to the spectra in each epoch, we generate a new random realization of the flux calibration error, as explained in Sec. 2.2, that we apply to both, the noisy input and the noiseless output spectrum. This is done to increase the robustness of the network against systematics in the flux calibration.  \nSimilarly, in each epoch, we further augment the data sets by applying a small, random offset (additive term) and scale (multiplicative term), both drawn from the normal distributions  (0 , 0 . 01) and  (1 , 1 . 01) , respectively.", '3.4 Loss function': 'Let us denote the simulated noiseless flux with 𝐹 𝜆 and the noisy spectrum with ̂ 𝐹 𝜆 , whereas the denoised spectrum with 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] where 𝐷 is the multivalued function describing the denoiser and 𝜆 indexes the spectral pixels. In the training loop, we minimize the root-mean-square error (RMSE) with the usual definition of  \n = √ √ √ √ 1 𝑁 𝑁 ∑ 𝜆 =1 ( 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] -𝐹 𝜆 ) 2 , (2)  \nwhere 𝑁 is the number of spectral pixels.  \nWe decided to use RMSE as the loss function for two reasons. The first one is its interpretability in the comparison of noisy and noiseless stellar spectra. In the ideal case, when the denoiser efficiently reduces the noise of the input, without introducing any artifacts, the RMS error measures the amount of residual noise. The second reason is its sensitivity to larger  \ndifferences between the actual output and the ideal one but it is easier to optimize using stochastic gradient descent than, for instance, the 𝐿 0 norm, which is only sensitive to the spectral pixel with the largest error.  \nWe expected that the largest errors would appear at the cores of the deepest absorption lines, as those are the features with the largest variance in the training sample. On the other hand, we aimed for a noise reduction that is more or less uniform over the entire wavelength coverage and also takes shallower absorption features into account with approximately equal weight. The mean absolute error (MAE) turned out to be too unsensitive at the strong spectral lines.', '3.5 Training': 'We optimized the weights of the network with the stochastic gradient descent (SGD) algorithm (Robbins & Monro, 1951) using a custom learning rate schedule. SGD is an iterative method that aims to minimize the loss function by updating the weights of a model based on the gradient of the loss with respect to the weights.  \nNotably, our implementation of SGD did not utilize Nesterov momentum, which is often used to accelerate convergence in relevant directions and dampen oscillations (Nesterov, 1983; Sutskever, Martens, Dahl, & Hinton, 2013). We made this choice to maintain simplicity and control of the learning process that we wanted to drive only with an adaptive learning rate scheduler.  \nDue to the large number of outputs and the sigmoid activation function before the network output the training experienced the problem of vanishing gradients. With a training  \nFIGURE3 Training and validation loss during a 7 , 500 epoch long training process, starting from the 1 , 000 th epoch. All loss values shown here are the epoch averages for each spectrum of the training (validation) data set. The irregularly fluctuating line of the validation loss sits on top of the smoother curve of the training loss, indicating that our model was neither under nor over-fitted. The value of the adaptive learning rate is also shown here in green. The learning rate is a series of step function according to how the adaptive learning rate scheduler decreases its value during the training process. Note: Cutting the first 1 , 000 epochs from the plot served the purpose of better visualization of the loss values. The initial epochs were characterized by an extremely steep loss decrease and including them would have made the plot less readable.  \n<!-- image -->  \nbatch size of 250 , the typical values of the gradient for the first layer happened to be as small as 10 -6 -10 -7 on average. To mitigate this problem without batch normalization, we trained our model with no L2 regularization (weight decay) and set the learning rate to an unusually high starting value of 100 . While this is generally discouraged due to the instabilities it can introduce into the loss function, we found that this unconventional technique greatly helped our model to converge.  \nDuring training, we adaptively reduced the learning rate when the model stopped improving for a predefined number of epochs, as measured from the validation loss. After some patience period, the learning rate was halved as seen in Fig. 3 . This approach allows the model to escape local minima and to fine-tune its parameters as the training progresses, ultimately leading to a more robust and generalized solution. In our case, the value of the patience parameter was set to 150 epochs.', '4 RESULTS': 'In this section we analyze the results obtained by training a denoising autoencoder neural network as described in Sec. 3.', '4.1 Measures of goodness': 'We measure the absolute and relative error after denoising by comparing the denoised spectrum 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] to the original simulated spectrum without any additional noise 𝐹 𝜆 , where ̂ 𝐹 𝜆 is the noisy spectrum. The definition of the absolute error for each pixel is  \nwhereas the relative error is taken to be  \n| 𝐸 𝜆 | = | | | 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] -𝐹 𝜆 | | | , (3)  \n| 𝑅 𝜆 | = | | | 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] -𝐹 𝜆 | | | 𝐹 𝜆 . (4)  \nWhen we characterize the denoising performace with a single number, we simply take the mean of | 𝐸 𝜆 | or | 𝑅 𝜆 | over every wavelength pixel.  \nTo provide information on the bias, we also report error metrics without taking the absolute value of the difference 𝐷 𝜆 [ ̂ 𝐹 𝜆 ] -𝐹 𝜆 . We denote these metrics without the absolute value as 𝐸 𝜆 and 𝑅 𝜆 .', '4.2 Denoising performance': 'As an example, in Fig. 4 and 5 , we plot the spectrum of a typical G giant and its reconstruction from a noisy simulated observation, as well as the residual and the relative error calculated as the ratio of the residual and the noiseless flux in each spectral pixel. As one would expect, the error is larger in the cores of the strong absorption lines, such as the lines of the calcium triplet because stellar spectra show a larger variance in the lines than in the continuum.  \nTo get a more general view, we plot the maximum and mean of the absolute errors | 𝐸 𝜆 | and relative errors | 𝑅 𝜆 | as functions of the signal-to-noise in Fig. 6 as well as functions of the fundamental stellar parameters in Fig. 7 -9 . We use the absolute values of the error metrics to generate data for these figures. This approach prevents the cancellation of positive and negative values in 𝐸 𝜆 and 𝑅 𝜆 , thereby preserving information on the magnitude of the errors.  \nUnsurprisingly, as visible in Fig. 6 , the performance of denoising scales with the median 𝑆 / 𝑁 of the spectra but the means of both the absolute error and the relative error are within a decade between 10 < 𝑆 / 𝑁 < 120 . The errors show a significant upturn toward lower 𝑆 / 𝑁 which limits the applicability of out network for really noisy spectra.  \nWhen looking at the residual noise as a function of 𝑇 eff in Fig. 7 , one can observed that the amount of remaining noise depends only slightly on temperature except for 𝑇 eff < 3750 K where the relative error shows a significantly larger spread than elsewhere. This can easily be attributed to the complexity of the optical spectra of M dwarfs due to the strong blanketing by molecular lines.  \nSome dependence of the residual error is visible as a function of log 𝑔 in Fig. 8 . The mean of the relative error is a few times larger toward the most extreme values of the validation set which is a typical result for networks trained and validated on data sets with the same parameter range coverage. The dependence of the residual noise on [Fe/H] , on the other hand, shows a more peculiar behaviour as it can be seen in Fig. 9 . Higher metallicity stellar spectra are reconstructed by the denoiser network with larger residual noise, very likely due to the presence of more and stronger absorption lines. Interestingly, the mean relative error scales more steeply with metellicity than the absolute error, hinting that a large number of spectral pixels is responsible for the increase in overall noise residual.', '4.3 Tests of generalizability': "To assess the robustness and generalizability of our denoising autoencoder, we took noiseless simulated spectra and added new, strong, artificial absorption lines with a Voigt profile and also increased the depth of some existing absorption lines. After adding the artificial lines, we added noise to the spectra and processed them with the denoising autoencoder that was trained with realistic absorption lines only.  \nOur quite plausible hypothesis was that the denoising autoencoder, more specifically its decoder, works by learning how to generate synthetic spectra from the latent space (i.e. global features) as opposed to learning how to remove the noise and reconstruct the continuum and spectral lines from local features.  \nIn Fig. 10 , we plot the calcium triplet region of a G giant star with artificially introduced additional absorption features. As the green curve of Fig. 10 indicates, the depth of the 8542 Å feature was increased whereas a completely new feature was introduced at 8600 Å. Fig. 10 shows the denoised spectrum in orange which exhibits the original line depth for the 8498 Åline and completely lacks the artificially introduced extra feature.  \nThis behaviour indicates that instead of simply performing pixel-wise denoising and returning the artificially introduced features without noise, the model learns the fundamental parameters from the simulated spectra and reconstructs the denoised spectrum from those, instead of 'blind' denoising of the input.", '5 CONCLUSION': 'We found that even simple neural networks, specifically denoising autoencoders (DAEs) consisting of fully connected  \nFIGURE 4 Top panel: Comparison of the original noiseless simulated stellar spectrum (solid blue line) and the denoised spectrum (dashed orange line) using our denoising autoencoder. Note that the two lines that represent the original and the denoised spectrum are almost entirely overlap, indicating high accuracy of our model. Bottom panel: The relative error calculated as the fraction of pixel-wise residual noise and the original noiseless flux. The mean and maximum of the relative error, as well as the Kullback-Leibler divergence of the original noiseless spectrum and the denoised spectrum are also indicated. Note, that the peaks of the relative error are aligned with the most prominent absorption features.  \n<!-- image -->  \nlayers with a bottleneck are capable of learning the characteristics of simulated stellar spectrum observations, We demonstrated that this capability enables the autoencoders to efficiently denoise simulated spectroscopic observation. For a wide range of fundamental stellar parameters, denoising can be performed at a sub-per-cent accuracy at 𝑆 / 𝑁 ≈ 10 .  \nBy introducing artificial absorption features into the simulated spectra, we were able to prove that the denoiser works by representing physical spectra in a low-dimensional latent space and reconstructs the output based on that instead of removing the noise by some kind of effective interpolation of local features. As a consequence, noise reduction capabilities depend on how well the encoder part can map the input noisy spectra into the latent space of the network bottleneck layer, as well as how accurately the decoder can transform the latent space variables into noiseless spectra. This knowledge can help training more complex networks with larger training sets that perform similarly well on noisier simulations and real observations.  \nFIGURE 5 A zoomed-in view of the calcium triplet for the same spectrum as seen in Fig. 4 . In addition to the relative error, the bottom panel also shows the pixel-wise absolute errors calculated as the absolute value of the difference between the denoised spectrum and the noiseless original spectrum.  \n<!-- image -->  \n<!-- image -->  \nS/N  \nFIGURE 6 Pixel-wise mean and maximum of the absolute and relative errors calculated from the residual noise of the denoised spectra for the entire validation data set as a function of median 𝑆 / 𝑁 . See the text for discussion.FIGURE 7 Pixel-wise mean and maximum of the absolute and relative errors calculated from the residual noise of the denoised spectra for the entire validation data set as a function of 𝑇 eff . See the text for discussion.  \n<!-- image -->  \n]  \n%  \n[  \n)  \n|  \nR  \n|  \n(  \ns  \nr  \no  \nr  \nr  \nE  \ne  \nv  \ni  \nt  \na  \nl  \ne  \nR  \nMaximum  \nMean  \n10  \n10  \n10  \n1  \n0  \n1  \n10  \n1  \n10  \n2  \n<!-- image -->  \nFIGURE 8 Pixel-wise mean and maximum of the absolute and relative errors calculated from the residual noise of the denoised spectra for the entire validation data set as a function of log 𝑔 . See the text for discussion.  \n<!-- image -->  \n<!-- image -->  \nFIGURE 9 Pixel-wise mean and maximum of the absolute and relative errors calculated from the residual noise of the denoised spectra for the entire validation data set as a function of [Fe/H] . See the text for discussion.  \n<!-- image -->', 'ACKNOWLEDGMENTS': 'This work is supported by the generosity of Eric and Wendy Schmidt, by recommendation of the Schmidt Futures program, the Ministry of Innovation and Technology NRDI Office grants OTKA NN 129148, the KDP-2021 program of the Ministry of Innovation and Technology from the source of the NRDI fund and by the the European Union project RRF-2.3.1-212022-00004 within the framework of the MILAB Artificial Intelligence National Laboratory.  \nThe simulation of all stellar spectrum observations used in this research, and the training of the DAE networks were completed on the Elephant and Volta servers at the Department of Physics & Astronomy of the Johns Hopkins University. Analysis was carried out on the Volta server, as well as on the Ampere A100 GPU server of the Wigner Scientific Computing Laboratory (WSCLAB) at the HUN-REN Wigner Research Centre for Physics, Hungary.', 'REFERENCES': 'Bohlin, R. C. (2014, June), AJ , 147 (6), 127.  \nFIGURE 10 Top panel: Comparison of a noiseless simulated spectrum (blue line) with artificially introduced extra absorption features (green line) to its denoised version (orange line). The denoised spectrum very precisely covers the original spectrum without the artificially added absorption lines. Bottom panel: The relative error of the denoising calculated as the pixel-wise absolute difference between the original noiseless spectrum with the artificial absorption features and the denoised spectrum, normalized by the noiseless flux. The relative error clearly indicates that the denoiser is unable to reconstruct the artificially introduced features which cause the peaks in relative error. Interestingly, additional peaks appear in the relative error at absorption lines which were not artificially modified, such as the two longer wavelength lines of the calcium triplet. This can be attributed to the strong correlation between the calcium lines that the denoising network learned from the training data and how the denoiser fails to benefit from the correlation when the depth of one of the lines is artificially changed.  \n<!-- image -->  \nBohlin, R. C., Mészáros, S., Fleming, S. W., Gordon, K. D., Koekemoer, A. M., & Kovács, J. (2017, May), AJ , 153 (5), 234.  \nCardelli, J. A., Clayton, G. C., & Mathis, J. S. (1989, October), ApJ , 345 , 245.  \nDESI Collaboration, Abareshi, B., Aguilar, J. et al. (2022, November), AJ , 164 (5), 207.  \nDey, A., Najita, J. R., Koposov, S. E. et al. (2023, February), ApJ , 944 (1), 1.  \nHinton, G. E., & Salakhutdinov, R. R. (2006), Science , 313 (5786), 504-507.  \nHirata, C. M., Gehrels, N., Kneib, J.-P., Kruk, J., Rhodes, J., Wang, Y., & Zoubian, J. (2012, April), arXiv e-prints , arXiv:1204.5151. Kramer, M. A. (1991), AIChE Journal , 37 (2), 233-243.  \nFIGURE 11 A zoomed-in view of the calcium triplet of the same spectra as in Fig. 10 . The bottom panel also includes a comparison of the pixel-wise relative and absolute errors of denoising.  \n<!-- image -->  \nKurucz, R. L. (1979, May), ApJS , 40 , 1-340.  \nMészáros, S., Allende Prieto, C., Edvardsson, B. et al. (2012, October), AJ , 144 (4), 120.  \nNesterov, Y. (1983), Doklady AN USSR , 269 , 543-547.  \nPaszke, A., Gross, S., Massa, F. et al. (2019), PyTorch: An Imperative Style, High-Performance Deep Learning Library. In H. Wallach, H. Larochelle, A. Beygelzimer, et al. (Eds.), Advances in Neural Information Processing Systems 32 pp. 8024-8035. Curran Associates, Inc.  \nRobbins, H., & Monro, S. (1951), The Annals of Mathematical Statistics , 22 (3), 400 - 407.  \nSugai, H., Tamura, N., Karoji, H. et al. (2014, July), Progress with the Prime Focus Spectrograph for the Subaru Telescope: a massively multiplexed optical and near-infrared fiber spectrograph. In S. K. Ramsay, I. S. McLean, & H. Takami (Eds.), Groundbased and Airborne Instrumentation for Astronomy V Vol. 9147, p. 91470T.  \nSutskever, I., Martens, J., Dahl, G., & Hinton, G. (2013, 17-19 Jun), On the importance of initialization and momentum in deep learning. In S. Dasgupta & D. McAllester (Eds.), Proceedings of the 30th International Conference on Machine Learning Vol. 28, pp. 1139-1147. Atlanta, Georgia, USA: PMLR.  \nTakada, M., Ellis, R. S., Chiba, M. et al. (2014, February), PASJ , 66 (1), R1.  \nTamura, N., Takato, N., Shimono, A. et al. (2016, August), Prime Focus Spectrograph (PFS) for the Subaru telescope: overview, recent progress, and future perspectives. In C. J. Evans, L. Simard, &H.Takami (Eds.), Ground-based and Airborne Instrumentation for Astronomy VI Vol. 9908, p. 99081M.  \nVincent, P., Larochelle, H., Bengio, Y., & Manzagol, P.-A. (2008), Extracting and composing robust features with denoising autoencoders. In International Conference on Machine Learning.'}