ABSTRACT

Several analyses of particle observations aim to determine the distribution functions of physical parameters that characterize observed systems. Some standard analysis methods determine these distributions by fitting mathematical models to the data. The accuracy of the fitting techniques depends on the treatment of the observations and their uncertainties. Here, we evaluate the performance of three fitting techniques by applying them to simulated electron observations, which are governed by the Poisson distribution. We specifically examine and quantify the accuracy of two standard chi-squared minimization techniques and a maximum-likelihood method. The chi-squared minimization techniques simplify the analysis by treating the measurement uncertainties as Gaussian errors. Although such a simplification reduces the complexity of the calculations in some occasions, it may lead to systematic errors in the determined parameters. On the other hand, the maximum-likelihood method considers the exact Poisson probability for each data-point and returns accurate parameters for all the examples we examine here. We highlight the importance of using the appropriate method when the observations are accompanied by significant statistical uncertainty. Nevertheless, the methods we examine here, converge to the same answer as the statistical uncertainty of the observations reduces.

1 INTRODUCTION

1.1 General background

Several analyses aim to estimate the parameters that characterize an observed system, by fitting a mathematical model to the observations. In such cases, the analytical model is a function of the system parameters and the parameters that most probably represent the observed system are those for which the model achieves the best fit. However, the accuracy of the returned parameters depends on the methodology that is used to determine the best fit to the observations. Typical fitting methods require the knowledge of the statistical uncertainty of the observations to optimize the fit, and thus, to determine the parameters (properties) of the observed system. Some methods adopt approximations regarding the statistical uncertainty of the observations and these approximations may lead to erroneous results. Previous studies have investigated the accuracy of standard fitting methods to measurements with statistical uncertainty governed by the Poisson distribution, such as Thomson scattering data (e.g. Stoneking & Den Hartog 1997), soft X-ray measurements (e.g. Yamada et al. 2019), and |$\gamma$|-ray spectra (e.g. Utsunomiya et al. 2018). The studies show that the methods that assume Gaussian statistical errors lead to erroneous results and the accuracy depends on the statistical uncertainty of the signal. In this study, we investigate the accuracy of similar fitting methods when used for the analysis of space particle observations by typical plasma instruments (e.g. McComas et al. 2004; Wilson 2015; Elliott et al. 2016). We specifically focus on applications that attempt to determine the plasma bulk parameters by optimizing a fit to the observations (e.g. Wilson, Bagenal & Persoon 2017; Scherer et al. 2022) and we evaluate the performance of each method as a function of the maximum number of the expected counts to be observed.

1.2 Statistical uncertainty of particle observations

We consider a detector that measures the number of particles that cross its detection area, within a specific time interval. For a steady particle flow, the probability |$P_{C}$| to get an individual measurement of C particles in the interval, is given by the Poisson distribution (e.g. Bevington & Robinson 2003):

(1)

where |$C_{\mathrm{exp}}$| is the expected (average) number of detected particles to be measured in the specific time interval. With some algebra, someone can show that the standard deviation of the Poisson distribution is |$\sqrt{C_{\mathrm{exp}}}$|⁠. However, we have to keep in mind that for small |$C_{\mathrm{exp}}$|⁠, the Poisson distribution is not symmetric and thus, the standard deviation cannot be used as a Gaussian uncertainty of the expected value. This study quantifies the systematic inaccuracies of the determined plasma bulk parameters, resulting by such simplifications regarding the statistical uncertainty of plasma particle observations.

2 FITTING METHODS

We examine three different fitting methods. The first two methods perform chi-square minimization, each using a different approach regarding the statistical uncertainty. The third method performs maximum-likelihood estimation, which optimizes the model by maximizing the likelihood of the observed data assuming they are characterized by a Poisson distribution.

2.1 Method A: |$\chi ^{2}$| minimization with data-determined uncertainty

Consider that we measure the number of particles as a function of their physical parameter X. In practice, we can sample the number of particles C in N discrete values of X. Then, our sampled distribution of counts is |$C(X_{i})$| with the subscript i denoting individual measurement samples. It is often convenient to assign a statistical uncertainty to each measurement sample |$\sigma \left[C(X_{i})\right]$|⁠. One approach is to consider |$\sigma \left[C(X_{i})\right]=\sqrt{C(X_{i})}$|⁠, as a direct outcome of the Poisson statistics. This is a convenient approach, as it assigns a statistical error to each individual measurement which is then easily propagated to any product of the measured counts (e.g. Nicolaou et al. 2015; Abraham et al. 2022). However, as discussed extensively in Nicolaou et al. (2020b), this approach considers that each sample |$C(X_{i})$| represents the expected value of the Poisson distribution at |$X_{i}$|⁠, but in reality, |$C(X_{i})$| is one event that may have resulted from a distribution of events with a different expected value. Nevertheless, when this approach is adopted, the metric to be minimized is

(2)

where |$M(\boldsymbol{p};X_{i})$| is the value of the model (model of the expected value) at each |$X_{i}$|⁠, and |$\boldsymbol{p}$| is the array of the model parameters to be optimized by the fitting. The users of this approach must be aware that |$C(X_{i})$| = 0 counts cannot be treated as data with zero uncertainty (⁠|$\sqrt{0}$|⁠), as this will result in |$\chi _{\mathrm{d-d}}^2=\infty$|⁠. Nicolaou et al. (2020b) show that in measurement samples with significant statistical uncertainty (determined by the peak of the counts and the samples of zero counts in a measurement), the accuracy of this method increases by including |$C(X_{i})$| = 0 counts in the analysis and assigning them an uncertainty |$\sigma \left[C(X_{i})\right]$| = 1 count, which is also suggested by Wilson (2015).

2.2 Method B: |$\chi ^{2}$| minimization with model-determined uncertainty

Another approach is to use the square root of the model values in each pixel as weights to the squared differences between the data and the model. The metric to be minimized is then (e.g. Tanabashi et al. 2018):

(3)

This approach considers a statistical uncertainty that is centred around the model values, which seems more appropriate for fitting measurements arriving from a distribution with expected value |$M(\boldsymbol{p};X_{i})$| and thus, a standard deviation |$\sqrt{M(\boldsymbol{p};X_{i})}$|⁠. However, this method still fails to handle non-Gaussian uncertainties and thus, is expected to lead to inaccuracies when the number of counts is small.

2.3 Method C: maximum-likelihood fitting

We also investigate the accuracy of another approach, introduced by Stoneking & Den Hartog (1997). The model is optimized by maximizing the probability of the observed counts being outcomes of the Poisson distribution with the expected value given by the model. According to equation (1), the probability of measuring |$C(X_{i})$| in each |$X_{i}$|⁠, while the expected value is |$M(\boldsymbol{p};X_{i})$|⁠, is given by

(4)

and thus, we can define the total probability for the fitted data set as

(5)

Therefore, in this approach, the model is optimized by searching for the set of parameters |$\boldsymbol{p}$| that maximize |$P_{t}$| or the logarithm of |$P_{t}$|⁠:

(6)

3 SIMULATIONS

We would like to investigate the accuracy of the three different fitting methods presented in Section 2, when applied to space plasma observations by typical space plasma instruments. We use the forward model method (e.g. Wilson et al. 2008, 2012; Cara et al. 2017; Criton, Nicolaou & Verscharen 2020) to simulate measurements by an electrostatic analyser instrument (e.g. Johnstone et al. 1997; McComas et al. 2017; Owen et al. 2020). We then apply the three fitting methods to the simulated observations and compare the results returned by each method to the corresponding parameters we use to simulate the observations. This comparison quantifies the accuracy of each method. Here, we model space plasma electrons with their energies following the isotropic kappa distribution function (e.g. Livadiotis & McComas 2011, 2013):

(7)

where |$n_{\mathrm{in}}$| is the plasma electron density, |$T_{\mathrm{in}}$| is the electron temperature, |$\kappa _{\mathrm{in}}$| is the kappa index, and m is the electron mass. |$E_{\mathrm{b,in}}$| is the bulk energy of the particles, |$\Omega$| is the angle between the individual particle velocity vector and the bulk velocity vector, and |$\Gamma$| is the gamma function. For simplicity, we consider that the bulk velocity vector of the plasma electrons is along |$\Theta =\Phi =0^{\circ }$|⁠. We further simplify by considering only 2D measurements, obtained at elevation angle |$\Theta =0^{\circ }$| (Nicolaou et al. 2020b). In this case, |$\Omega$| is given by the azimuth angle |$\Phi$| and the energy (velocity) distribution function becomes:

(8)

For our study, we simulate measurements of such electron plasma by an electrostatic analyser which measures the number of particles in 64 discrete energies E ranging from 1 eV to 5 keV (or speeds |$V=\sqrt{\frac{2E}{m}}$| ranging from |$\sim$|600 km|$\, \mathrm{s}^{-1}$| to |$\sim$|42 000 km|$\, \mathrm{s}^{-1}$|⁠), and 32 azimuth |$\Phi$| directions ranging from 0|$^{\circ }$| to 360|$^{\circ }$|⁠. This design is based on the Solar Wind Analyser–Electron Analyser System on board Solar Orbiter (Owen et al. 2020). The expected number of counts in each E, |$\Phi$| pixel is approximately (e.g. Nicolaou, Livadiotis & Wicks 2020c; Nicolaou 2023):

(9)

where |$G_{\mathrm{F}}$| is the effective geometric factor of the instrument (for detailed definition see Collinson et al. 2012) and |$\Delta \tau$| is the acquisition time for each individual measurement. For our simulations, we set |$\Delta \tau$| = 1 ms. The simulation takes into account that the number of counts registered in each pixel follows the Poisson distribution with expected (mean) value |$C_{\mathrm{exp}}(E,\Theta =0^{\circ },\Phi )$|⁠, governed by the probability given in equation (1) (e.g. Nicolaou et al. 2020a; Zhang, Verscharen & Nicolaou 2024).

We evaluate the fitting methods we describe in Section 2 by applying them to simulated observations with |$n_{\mathrm{in}}$| =10 cm−3, |$E_{\mathrm{b,in}} \sim$| 0.71 eV, |$k_{\mathrm{B}}T_{\mathrm{in}}$| = 12 eV, and |$\kappa _{\mathrm{in}}$| = 3. However, to examine the accuracy of each method as a function of the measurement uncertainty, we simulate measurements of the same plasma but for different values of |$G_{\mathrm{F}}$|⁠. As it can be seen by equation (9), a change in |$G_{\mathrm{F}}$| results in an equivalent change in |$C_{\mathrm{exp}}$|⁠. Here, we run simulations for 18 different values of |$G_{\mathrm{F}}$|⁠, spanning from |$\sim 8.1\times 10^{-9}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$| to |$\sim 8.1\times 10^{-8}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$|⁠, which correspond to a maximum |$C_{\mathrm{exp}}$| ranging from 11.3 to 113 counts. For each |$G_{\mathrm{F}}$| value, we simulate 1000 observations. We note that the 1000 simulated samples for the same plasma and |$G_{\mathrm{F}}$|⁠, are different from each other since our simulations take into account that the statistical variations follow the Poisson distribution. Each observation sample is analysed by each of the three methods, such that, each fitting determines the plasma density |$n_{\mathrm{out}}$|⁠, temperature |$T_{\mathrm{out}}$|⁠, and kappa index |$\kappa _{\mathrm{out}}$|⁠. Note that the plasma bulk energy is not a free parameter here, considering that in typical solar wind applications it can be determined precisely from proton observations. Unlike electrons, proton energy distribution functions have a distinct peak at the bulk energy, due to their small thermal-to-bulk energy ratio (see discussions by Nicolaou et al. 2020b, 2021). Then, for each |$G_{\mathrm{F}}$| and each method, we determine 1000 |$n_{\mathrm{out}}$|⁠, temperature |$T_{\mathrm{out}}$|⁠, and kappa index |$\kappa _{\mathrm{out}}$| values and we examine their histograms, median values, lower, and higher quartiles.

The middle panel of Fig. 1 shows the 2D simulated observations for |$G_{\mathrm{F}}\sim 8.1\times 10^{-9}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$|⁠, which is the lowest |$G_{\mathrm{F}}$| value we examine in this paper. The maximum number of |$C_{\mathrm{exp}}$| for the specific case is |$\sim$|11.3 counts and the measurement is characterized by significant statistical uncertainty, which is evident from the variation of the counts resulting in a ‘blurred’ distribution of counts over E and |$\Phi$|⁠. The left and right panels show 1D cuts of C as a function of E for the two azimuth sectors around |$\Phi = 180^{\circ }$|⁠, and the two azimuth sectors around |$\Phi = 0^{\circ }$|⁠, respectively. On the same panels, we show the corresponding 1D cuts of the models optimized by each of the three different methods we investigate. As it is theoretically expected (Livadiotis 2007), each method clearly optimizes the fit model differently and thus, by determining different plasma bulk parameters. In the next section, we present the returned parameters by each method when applied to our simulations and quantify the accuracy of each method as a function of the highest number of expected counts.

Simulated observations and fit models. (Middle) 2D array of simulated counts of electrons with $n_{\mathrm{in}}$ =10 cm−3, $E_{\mathrm{b,in}} \sim$ 0.71 eV, $k_{\mathrm{B}}T_{\mathrm{in}}$ = 12 eV, $\kappa _{\mathrm{in}}$ = 3, and an instrument with $G_{\mathrm{F}}\sim 8.1\times 10^{-9}$  $\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$. (Left) 1D cuts of the modelled observations (black circles) and the three model fits (coloured lines) for two azimuth sectors around $\Phi = 180^{\circ }$, and (right) for the two azimuth sectors around $\Phi = 0^{\circ }$.
Figure 1.

Simulated observations and fit models. (Middle) 2D array of simulated counts of electrons with |$n_{\mathrm{in}}$| =10 cm−3, |$E_{\mathrm{b,in}} \sim$| 0.71 eV, |$k_{\mathrm{B}}T_{\mathrm{in}}$| = 12 eV, |$\kappa _{\mathrm{in}}$| = 3, and an instrument with |$G_{\mathrm{F}}\sim 8.1\times 10^{-9}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$|⁠. (Left) 1D cuts of the modelled observations (black circles) and the three model fits (coloured lines) for two azimuth sectors around |$\Phi = 180^{\circ }$|⁠, and (right) for the two azimuth sectors around |$\Phi = 0^{\circ }$|⁠.

4 RESULTS

Fig. 2 shows histograms of the parameters determined by each of the three methods when applied to 1000 samples of simulated measurements with instrument |$G_{\mathrm{F}}\sim 8.1\times 10^{-9}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$| (maximum |$C_{\mathrm{exp}}$| = 11.3). The left panel shows histograms of the plasma density |$n_{\mathrm{out}}$|⁠, the middle panel shows histograms of |$k_{\mathrm{B}}T_{\mathrm{out}}$|⁠, while the right panel shows histograms of |$\kappa _{\mathrm{out}}$|⁠. It is clear that the histograms of the plasma parameters determined by the maximum-likelihood method (Method C, blue) are more representative of the actual plasma parameters (input parameters shown by the vertical magenta), compared to the parameters determined by the chi-squared minimization with data-driven uncertainty (Method A, grey) and the chi-squared minimization with model-driven uncertainty (Method B, red). The |$n_{\mathrm{out}}$| determined by Method A has a median value |$\sim$| 7.7 cm−3, which is more than 20 per cent smaller than the corresponding input density. The median of |$n_{\mathrm{out}}$| determined by Method B is |$\sim$| 11.9 cm−3, almost 20 per cent larger than the input value. Importantly, Method A and Method B fail to recover the actual density within the range of their |$n_{\mathrm{out}}$| histograms. On the other hand, the median density by Method C recovers the actual density (⁠|$n_{\mathrm{in}}$|⁠) within 2 per cent. The median value of |$k_{\mathrm{B}}T_{\mathrm{out}}$| determined by Method A is |$\sim$|11.5 eV and its higher quartile |$\sim$|11.7 eV which are both lower than |$k_{\mathrm{B}}T_{\mathrm{in}}$| = 12 eV. Method B overestimates the plasma temperature by |$\sim$|12.5 per cent since the corresponding histogram has a median value |$\sim$|13.8 eV and a lower quartile |$\sim$|13.5 eV. Method C is recovering almost perfectly the input temperature since it returns a median value of |$\sim$|12 eV, a lower quartile |$\sim$|11.8 eV, and a higher quartile |$\sim$|12.2 eV. Method A slightly overestimates the kappa index since its |$\kappa _{\mathrm{out}}$| histogram has a median |$\sim$|3.2. The left tail of the histogram, however, includes the input value |$\kappa _{\mathrm{in}}$| = 3. Method B underestimates the kappa index as the returned values have a median |$\sim$|2.3, and their range does not include the input value. The median value of |$\kappa _{\mathrm{out}}$| determined by Method C is almost identical to the input value, while the lower and higher quartiles are 2.95 and 3.06, respectively.

Histograms of plasma bulk parameters determined by applying the three different fit methods to the same measurement samples. (Left) Density $n_{\mathrm{out}}$ (middle) temperature $k_{\mathrm{B}} T_{\mathrm{out}}$, and (right) kappa index $\kappa _{\mathrm{out}}$, determined by using Method A: chi-squared minimization with data-driven uncertainty (grey), Method B: chi-squared minimization with model-driven uncertainty (red), and Method C: maximum-likelihood method (blue). The simulations consider $n_{\mathrm{in}}$ = 10 cm−3, $k_{\mathrm{B}} T_{\mathrm{in}}$ = 12 eV, $\kappa _{\mathrm{in}}$ = 3, and an instrument $G_{\mathrm{F}}\sim 8.1\times 10^{-9}$  $\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$.
Figure 2.

Histograms of plasma bulk parameters determined by applying the three different fit methods to the same measurement samples. (Left) Density |$n_{\mathrm{out}}$| (middle) temperature |$k_{\mathrm{B}} T_{\mathrm{out}}$|⁠, and (right) kappa index |$\kappa _{\mathrm{out}}$|⁠, determined by using Method A: chi-squared minimization with data-driven uncertainty (grey), Method B: chi-squared minimization with model-driven uncertainty (red), and Method C: maximum-likelihood method (blue). The simulations consider |$n_{\mathrm{in}}$| = 10 cm−3, |$k_{\mathrm{B}} T_{\mathrm{in}}$| = 12 eV, |$\kappa _{\mathrm{in}}$| = 3, and an instrument |$G_{\mathrm{F}}\sim 8.1\times 10^{-9}$|  |$\mathrm{m}^2\, \mathrm{sr}\, \mathrm{eV}\, \mathrm{eV}^{-1}$|⁠.

Fig. 3 shows the median values of |$n_{\mathrm{out}}$| (top), |$k_{\mathrm{B}} T_{\mathrm{out}}$| (middle), and |$\kappa _{\mathrm{out}}$| (bottom), determined by each method when applied to simulations of the same plasma (⁠|$n_{\mathrm{in}}$| = 10 cm−3, |$E_{\mathrm{b,in}} \sim$| 0.71 eV, |$k_{\mathrm{B}}T_{\mathrm{in}}$| = 12 eV, |$\kappa _{\mathrm{in}}$| = 3), but for 18 different values of |$G_{\mathrm{F}}$|⁠, and thus, 18 different values of maximum |$C_{\mathrm{exp}}$|⁠. Each data point is the median value and the shaded areas indicate the range bounded by the lower and higher quartiles of the determined parameters (from the analysis of 1000 samples for each |$G_{\mathrm{F}}$| value). The results show that Method C determines all the plasma parameters with significant accuracy for all the |$G_{\mathrm{F}}$| values we examine. We also see that all three methods converge to the right answer as |$C_{\mathrm{exp}}$| increases. Additionally, for the entire range of |$G_{\mathrm{F}}$| range we examine, Method A underestimates the plasma density while Method B overestimates it. Method B, however, seems to provide more accurate |$n_{\mathrm{out}}$| values for the entire range of |$G_{\mathrm{F}}$| (top panel of Fig. 3). On the other hand, the plasma temperature determined by Method A, is more accurate than that determined by Method B. The median of |$k_{\mathrm{B}} T_{\mathrm{out}}$| determined by Method A is always within 1 eV of the actual value, while the corresponding median determined by Method B is larger than the actual value by more than 1 eV for all the cases with maximum |$C_{\mathrm{exp}}$| < 35 counts. The kappa index is overestimated by Method A. The median of |$\kappa _{\mathrm{out}}$| determined by Method A approaches its peak |$\sim$|3.5 (0.5 larger than |$\kappa _{\mathrm{in}}$|⁠) at maximum |$C_{\mathrm{exp}} \sim$| 35 counts. Method B, on the other hand, underestimates |$\kappa _{\mathrm{out}}$| by more than 0.5 for the cases with maximum |$C_{\mathrm{exp}}$| < 30 counts. Finally, we see that the medians of |$k_{\mathrm{B}} T_{\mathrm{out}}$| and |$\kappa _{\mathrm{out}}$| determined by Method A, do not vary monotonically with maximum |$C_{\mathrm{exp}}$|⁠. In the next section, we discuss these results and their implications.

Plasma bulk parameters determined by each of the examined methods as a function of the maximum number of expected counts. (Top) $n_{\mathrm{out}}$, (middle) $k_{\mathrm{B}} T_{\mathrm{out}}$, and (bottom) $\kappa _{\mathrm{out}}$ versus maximum $C_{\mathrm{exp}}$, as determined by Method A (grey), Method B (red), and Method C (blue).
Figure 3.

Plasma bulk parameters determined by each of the examined methods as a function of the maximum number of expected counts. (Top) |$n_{\mathrm{out}}$|⁠, (middle) |$k_{\mathrm{B}} T_{\mathrm{out}}$|⁠, and (bottom) |$\kappa _{\mathrm{out}}$| versus maximum |$C_{\mathrm{exp}}$|⁠, as determined by Method A (grey), Method B (red), and Method C (blue).

5 DISCUSSION AND CONCLUSIONS

We have examined the accuracy of three standard fitting methods when applied to space plasma observations by a typical electrostatic analyser. We evaluate the accuracy of the methods by analysing simulated plasma electron observations. We focus mostly on quantifying the accuracy of these methods when applied to observations with significant statistical uncertainty (maximum |$C_{\mathrm{exp}}$| < 150 counts).

We highlight that the maximum-likelihood fitting (Method C) returns all the plasma parameters with remarkable accuracy, even when applied to the samples with the largest uncertainty we consider here. The superiority of the maximum-likelihood fitting over chi-squared minimization fitting methods was also shown by other studies considering different physical systems (e.g. Stoneking & Den Hartog 1997; Utsunomiya et al. 2018; Yamada et al. 2019). Our study is specifically quantifying the results for space plasma measurements. The example with the smallest peak we examine here has maximum |$C_{\mathrm{exp}} \sim$| 11.3 counts and thus, a relative uncertainty (standard deviation over mean value) |$1/\sqrt{11.3} \sim$| 30 per cent. It is then expected that the majority of the data points have uncertainty larger than 30 per cent, since their values are below the peak value. Fig. 1 shows the specific modelled example. The large variation of the modelled counts reflects their significant statistical uncertainty. Nevertheless, the maximum-likelihood fitting is returning accurately the plasma bulk parameters and with a remarkably small statistical error (small range of returned parameters, see Figs 2 and 3).

We show that the classic  |$\chi ^2$| minimization methods lead to non-negligible, systematic inaccuracies. The systematic errors resulted from Method A, which uses |$\sqrt{C}$| as the statistical uncertainty of each data point, are consistent with the findings by Nicolaou et al. (2020b) who evaluate the specific method for a similar space plasma application. Their study also shows that there is an inter-dependence of the kappa index and the temperatures determined by the specific method, which is reflected in the apparent anticorrelation of the median values of |$k_{\mathrm{B}} T_{\mathrm{out}}$| and |$\kappa _{\mathrm{out}}$| in the middle and bottom panel of Fig. 3, respectively. Here, we conclude that although it is convenient, and occasionally not problematic, to consider |$\sqrt{C}$| as the statistical uncertainty of the observations (e.g. Nicolaou & Livadiotis 2020), it may lead to significant errors of fitting applications, especially when the overall number of counts in the analysed sample is small.

We also show that Method B, which weights the fitting by the square root of the model, determines more accurate densities than those determined by Method A. In Fig. 1, we see that the model optimized by Method B has higher values than the model optimized by Method A. This is expected because when using the |$\sqrt{C}$| weighting (Method A), the model ‘passes’ closer to data with smaller C, which are assigned a smaller uncertainty (⁠|$\sqrt{C}$|⁠). On the other hand, the |$\sqrt{M}$| weighting results to a model that aims to minimize its difference from smaller and higher values, regardless. This approach seems more appropriate for the core of the distribution, when the observed number of counts is larger and their statistical error is closer to a Gaussian error. This is possibly why, the densities returned by Method B are more accurate, considering that the density of the plasma is mostly determined by the accurate characterization of the distribution’s core.

Nevertheless, Method B is less accurate than Method A when it comes to the determination of the plasma temperature and kappa index. In Fig. 1, we show that the model by Method B has larger number of counts than any other model in the high energy regime, which possibly explains why Method B always returns higher |$k_{\mathrm{B}} T_{\mathrm{out}}$| and smaller |$\kappa _{\mathrm{out}}$| than the other two methods. The |$C_{\mathrm{exp}}$| at the tails of the counts distribution is considerably low, and thus, the Gaussian uncertainty approach fails completely. So, both Methods A and B are expected to lead to inaccuracies when characterizing the tails. However, in the very low |$C_{\mathrm{exp}}$| range, Method B will force a model to pass closer to the higher counts of the sample for which the weight is more relaxed. For example, consider an energy regime with |$C_{\mathrm{exp}}\ll$| 1 count. The measurement is then expected to contain samples with C = 0 counts, a fewer samples with C = 1 count, and even a fewer samples with higher number of counts. However, each chi-squared term by Method B is weighted by the square root of the model (equation 3). Therefore, the optimization will ‘force’ the model to pass closer to the higher counts in the sample even though they are much fewer, since a chi-squared term in equation (3) with non-zero |$C(X_{i})$| and extremely small |$M(X_{i})$| will be very large.

One advantage of Methods A and B compared to Method C, is their ability to perform calculations with a large number of counts. Typical computers have no problem to perform calculations of |$\chi _{\mathrm{d-d}}^2$| and |$\chi _{\mathrm{m-d}}^2$|⁠, even for C and |$C_{\mathrm{exp}}$| of several millions of counts. However, the Poisson probability calculated in Method C contains terms which may cause overflow errors for measurements of a few hundreds of counts. For instance, for |$C(X_{i}) = M(\boldsymbol{p};X_{i}) = 140$| counts, the |$M(\boldsymbol{p};X_{i})^{C(X_{i})}$| term in equation (4) is larger than |$10^{300}$|⁠. Therefore, for practical reasons, Methods A and B may be preferred over Method C for samples with large number of counts, for which all three methods converge to the same answer anyway.

Finally, we would like to highlight that in this study we use only one set of plasma parameters for our simulations. The performance of each method is expected to vary with the plasma parameters, but in any case, it can be evaluated by using the techniques we use in this paper. We also neglect the potential effects of noise in the simulated measurements, which have been investigated in previous studies (e.g. Nicolaou 2023; Wang et al. 2024).

ACKNOWLEDGEMENTS

GN thanks David J. McComas, Jesse Coburn, and Xiangyu Wu for the helpful discussions. GL was funded in part by the IBEX mission as part of NASA’s Explorer Program (80NSSC18K0237) and IMAP mission as a part of NASA’s Solar Terrestrial Probes (STP) Program (80GSFC19C0027). CI acknowledges STFC’s support through the PhD degree (ST/X508858/1).

DATA AVAILABILITY

This study uses simulated data which will be provided by the author upon request.

REFERENCES

Abraham
 
J. B.
 et al. ,
2022
,
ApJ
,
931
,
118
 

Bevington
 
P. R.
,
Robinson
 
D. K.
,
2003
,
Data Reduction and Error Analysis for the Physical Sciences
.
McGraw-Hill
,
New York

Cara
 
A.
 et al. ,
2017
,
J. Geophys. Res.: Space Phys.
,
122
,
1439
 

Collinson
 
G. A.
 et al. ,
2012
,
Rev. Sci. Instrum.
,
83
,
033303
 

Criton
 
B.
,
Nicolaou
 
G.
,
Verscharen
 
D.
,
2020
,
Appl. Sci.
,
10
,
8483
 

Elliott
 
H. A.
,
McComas
 
D. J.
,
Valek
 
P.
,
Nicolaou
 
G.
,
Weidner
 
S.
,
Livadiotis
 
G.
,
2016
,
ApJS
,
223
,
19
 

Johnstone
 
A. D.
 et al. ,
1997
, in
Escoubet
 
C. P.
,
Russell
 
C. T.
,
Schmidt
 
R.
, eds,
Peace: A Plasma Electron and Current Experiment
.
Springer
,
Dordrecht, the Netherlands
, p.
351
 

Livadiotis
 
G.
,
2007
,
Phys. A: Stat. Mech. Appl.
,
375
,
518
 

Livadiotis
 
G.
,
McComas
 
D. J.
,
2011
,
ApJ
,
741
,
88
 

Livadiotis
 
G.
,
McComas
 
D. J.
,
2013
,
Space Sci. Rev.
,
175
,
183
 

McComas
 
D. J.
 et al. ,
2004
,
J. Geophys. Res.: Space Phys.
,
109
,
A02104
 

McComas
 
D. J.
 et al. ,
2017
,
Space Sci. Rev.
,
213
,
547
 

Nicolaou
 
G.
,
2023
,
Astrophys. Space Sci.
,
368
,
3
 

Nicolaou
 
G.
,
Livadiotis
 
G.
,
2020
,
Entropy
,
22
,
541
 

Nicolaou
 
G.
,
McComas
 
D.
,
Bagenal
 
F.
,
Elliott
 
H.
,
Ebert
 
R.
,
2015
,
Planet. Space Sci.
,
111
,
116
 

Nicolaou
 
G.
,
Wicks
 
R. T.
,
Rae
 
I. J.
,
Kataria
 
D. O.
,
2020a
,
Space Weather
,
18
,
e2020SW002559
 

Nicolaou
 
G.
,
Wicks
 
R.
,
Livadiotis
 
G.
,
Verscharen
 
D.
,
Owen
 
C.
,
Kataria
 
D.
,
2020b
,
Entropy
,
22
,
103
 

Nicolaou
 
G.
,
Livadiotis
 
G.
,
Wicks
 
R. T.
,
2020c
,
Entropy
,
22
,
212
 

Nicolaou
 
G.
 et al. ,
2021
,
AA
,
656
,
A10
 

Owen
 
C. J.
 et al. ,
2020
,
A&A
,
642
,
A16
 

Scherer
 
K.
,
Husidic
 
E.
,
Lazar
 
M.
,
Fichtner
 
H.
,
2022
,
A&A
,
663
,
A67
 

Stoneking
 
M. R.
,
Den Hartog
 
D. J.
,
1997
,
Rev. Sci. Instrum.
,
68
,
914
 

Tanabashi
 
M.
 et al. ,
2018
,
Phys. Rev. D
,
98
,
030001
 

Utsunomiya
 
H.
 et al. ,
2018
,
Nucl. Instrum. Methods. Phys. Res. A
,
896
,
103
 

Wang
 
Z.
 et al. ,
2024
,
Adv. Space Res.
,
74
,
5282
 

Wilson
 
R. J.
,
2015
,
Earth Space Sci.
,
2
,
201
 

Wilson
 
R. J.
,
Tokar
 
R. L.
,
Henderson
 
M. G.
,
Hill
 
T. W.
,
Thomsen
 
M. F.
,
Pontius Jr
 
D. H.
,
2008
,
J. Geophys. Res.: Space Phys.
,
113
,
A12218
 

Wilson
 
R. J.
,
Delamere
 
P. A.
,
Bagenal
 
F.
,
Masters
 
A.
,
2012
,
J. Geophys. Res.: Space Phys.
,
117
,
A03212
 

Wilson
 
R. J.
,
Bagenal
 
F.
,
Persoon
 
A. M.
,
2017
,
J. Geophys. Res.: Space Phys.
,
122
,
7256
 

Yamada
 
S.
 et al. ,
2019
,
PASJ
,
71
,
75
 

Zhang
 
H.
,
Verscharen
 
D.
,
Nicolaou
 
G.
,
2024
,
Space Weather
,
22
,
e2023SW003671
 

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.