Flexible Two-point Selection Approach for Characteristic Function-based Parameter Estimation of Stable Laws

Stable distribution is one of the attractive models that well describes fat-tail behaviors and scaling phenomena in various scientific fields. The approach based upon the method of moments yields a simple procedure for estimating stable law parameters with the requirement of using momental points for the characteristic function, but the selection of points is only poorly explained and has not been elaborated. We propose a new characteristic function-based approach by introducing a technique of selecting plausible points, which could bring the method of moments available for practical use. Our method outperforms other state-of-art methods that exhibit a closed-form expression of all four parameters of stable laws. Finally, the applicability of the method is illustrated by using several data of financial assets. Numerical results reveal that our approach is advantageous when modeling empirical data with stable distributions.

Stable distribution is one of the attractive models that well describes fat-tail behaviors and scaling phenomena in various scientific fields. The approach based upon the method of moments yields a simple procedure for estimating stable law parameters with the requirement of using momental points for the characteristic function, but the selection of points is only poorly explained and has not been elaborated. We propose a new characteristic function-based approach by introducing a technique of selecting plausible points, which could bring the method of moments available for practical use. Our method outperforms other state-of-art methods that exhibit a closed-form expression of all four parameters of stable laws. Finally, the applicability of the method is illustrated by using several data of financial assets. Numerical results reveal that our approach is advantageous when modeling empirical data with stable distributions.

I. INTRODUCTION
A fundamental theory of stochastic processes in various scientific fields is the generalized central limit theorem (GCLT), which points out that the sum of independent and identically distributed random variables converges only to the family of stable distribution [1]. There are some challenges to overcome the analytic difficulties of stable distributions since the probability density function (PDF) is not always expressed in a closed form. Numerically approximated expressions are known in symmetric cases (β = 0) based on hypergeometric functions, but those in unrestricted asymmetric cases are often too complex for estimating the parameters of the stable distribution [2]. More practically, the estimation of all parameters is the most basic and necessary process for any application, but it remains to be one of the most controversial issues when attempting to detect stable laws. Numerous approaches have been studied for the parameter estimation. The primary approaches include the approximate maximum likelihood estimation [3][4][5][6], the bayesian based method [7], the quantile method (QM) [8,9], the fractional lower order moment (FLOM) method [10,11], the method of log-cumulant (MLOC) [12,13], the characteristic function-based (CF-based) method [14][15][16][17][18], and their hybrid combinations. Many of them tend to have different kinds of drawbacks, such as restrictions of parameter ranges, complex estimation algorithms, high computational costs, requirements of larger datasets, and low accuracy. To the best of our knowledge, the FLOM, MOLC, and QM and some class of the CF-based methods [15][16][17][18] provide closed-form estimators of stable laws.
The CF-based method is perhaps the largest classification group, including a variety of methods and approaches developed under different techniques. In particluar, Press (1972) [15] presents the method of moments, which offers a simple approach to estimate all four parameters of stable distribution using the characteristic function (CF) evaluated at four arbitrary points. The biggest advantage of this method is that it is likely to have less drawbacks compared to other primary methods, but it carries a fundamental problem. Without appropriate points given, the performance is poor, and unfortunately Press leaves unsolved the crucial idea about the choice of points at which the CF should be evaluated. The selection of the points has long been an open question, although several studies have made an effort to improve the method of moments by reducing the use of points from four to two and discussing their choice. Krutto (2016Krutto ( , 2018 provides some guidance on how the two positive points should be chosen through empirical searches relying on the cumulant function [17,18]. Bibalan et al. (2017) focus on the absolute value of the CF and suggest an algorithmic approach where a positive point is fixed for each scaling parameter [16]. They show accurate estimates within certain parameter ranges, but their method fails to support a wider range of parameter spaces. Thus, these approaches are not comprehensive, so that the method of detecting more appropriate points related to the CF is required for practical uses.
In this paper, we propose an effective and practical method for estimating stable laws. We greatly improve the method of moments by introducing a new technique for the selection of two positive points at which the CF is evaluated. The technique is developed over the extension of both algorithmic and empirical search approaches. The idea of empirical search plays a role in determining the scaling related estimates, which take crucial responsibility for indicating statistical values derived in the estimation process, whereas the concept of the algorithmic approach yields various ideas of inferences based on the absolute value of the CF. Our approach realizes the possibility of choosing different values of points depending on the index parameter α, which is a new perspective. We assess and compare the performance of our method to those of other methods in terms of the Mean Squared Error (MSE) criterion and the Kolmogorov-Smirnov (KS) distance. Our proposed method generally outperforms all the other state-of-art methods that exhibit closedform expressions for all four parameters of stable laws.

arXiv:2005.11499v1 [stat.ME] 23 May 2020
It is practically straightforward and assures that there is no restriction of parameter ranges, except for α = 1 due to the discontinuous form of the one-parameterization CF. Finally, we apply our method to price fluctuation behaviors of several financial assets to examine the appropriateness for practical uses.
This paper is organized as follows. Section 2 shows preliminaries on stable distribution and its basic properties. We follow in the next section to describe the existing methods for estimating the parameters of stable laws. In section 4 we propose a new technique of the CF-based parameter estimation method. The arguments for the selection of points at which the CF should be evaluated are discussed. In section 5 we report the performance with the comparison to other representative methods and present that our method provides accurate estimates of stable distribution. The last section shows application to financial data and confirms that our method is applicable for empirical studies.

II. STABLE DISTRIBUTION
In this section, we summarize the basis and properties of the stable distribution. We explain the definition of the stable distribution and its properties.

A. Basis of stable distribution
Stable distribution, also known as α-stable distribution, or Lévy's stable distribution, was first introduced by Paul Lévy (1937) [19], which is a family of parametric distribution with tails that are expressed as powerfunctions. In the far tails the PDF can be written as [20], and the cumulative distribution function (CDF) written as, where c α is a constant value [sin(πα/2)Γ(α)]/π. Stable distribution is represented by four parameters; the tail index parameter α ∈ (0, 2] representing the fatness of the tail, the skewness parameter β ∈ [−1, 1], the scaling parameter γ > 0, and the location parameter δ ∈ R. Especially the parameters α and β determine the shape of distribution, including various forms of widely-known distributions such as the Gaussian and Cauchy distribution. Smaller value of α indicates fatter tails and hence it is well known that the variance diverges for 0 < α < 2, and also the mean cannot be defined for 0 < α ≤ 1. Note that if β = 0, the distribution is symmetric, if β > 0, right-tailed, and if β < 0, left-tailed.
The definition of stable distribution is that the linear combination of independent random variables that follow a stable distribution with tail index α invariably becomes again a stable distribution with the same tail index α. More particularly, when variables X 1 , X 2 are i.i.d. copies of a random variable X and a, b are positive constant numbers, X is said to be stable and follows a stable distribution if there is a positive constant number c and a real number d ∈ R that satisfies also known for stability property. When a variable X follows a stable distribution, the notation X d = S(α, β, γ, δ) is often used, where d = denotes equality in distribution [21]. Variable X can be standardized according to the following property: Another important property of stable distribution is the GCLT, which implies that the only possible limit distributions for sums of i.i.d random variables is a family of stable distribution. When α = 2, that is, when i.i.d. random variables have finite variance, the limit distribution then becomes a Gaussian according to the well-known classical Central Limit Theorem (CLT).

B. Characteristic function
The PDF of stable distribution cannot be written in a closed form except for some cases; Cauchy distribution (α = 1, β = 0), Lévy distribution (α = 1/2, β = 1), and Gaussian distribution (α = 2). Alternatively, the features are expressed by the characteristic function (CF), ϕ(k), which is the Fourier transform of the PDF. By taking the inverse Fourier transform of the CF, the PDF can be obtained as When variable X follows a stable distribution with S(α, β, γ, δ), the CF is shown as which corresponds to the one-parameterization form of S(α, β, γ, δ; 1) in Nolan (2003) [22]. This is the most popular parameterization among many other forms of the stable distribution owing to the simplicity of the form. Figure 1 shows the standardized stable distributions with the one-parameterization form for different parameters of α and β, as an example. One-parameterization is preferred when one is interested in the basic properties of the distribution, but the CF takes a discontinuous form at α = 1. Nolan suggests the use of the zero-parameterization form S(α, β, γ, δ 0 ; 0) with different ω(k, α) shown as giving a more complex form, but provides a continuous form. The only difference between the parameterization is the location parameter, which they are related by In this paper, we employ the simple oneparameterization, as we are interested in estimating the four parameters through the CF, and many existing estimation methods comply with that form. However, since this CF does not have a continuous form at α = 1, arguments with different parameterizations may be more appropriate for discussing distributions when we already know that α is 1, for instance, the case of Cauchy distribution (α = 1, β = 0).

III. PARAMETER ESTIMATION OF STABLE LAWS
This section gives an overview of the methods for the parameter estimation of the stable distribution. We review two major methods, both of which are considered as an analytical approach that provides a closedform expression of the estimates-the quantile method and the characteristic function-based method (CF-based method). Several different approaches are explained for the CF-based method.

A. Quantile method
McCulloch (1986) proposes the use of five sample quantiles x 0.05 , x 0.25 , x 0.5 , x 0.75 , and x 0.95 as an informative measure for estimating the four parameters of stable laws, known as the quantile method (QM) [9]. He improves the former method of Fama & Roll (1971) by eliminating bias in estimates and relaxing estimation restrictions [8]. The idea is to calculate the functions φ i (α, β) (i = 1, 2, 3, 4), where the relationships between the function values and the parameters are already studied and known beforehand. The method first sets out to estimate α and β by using the functions φ 1 (α, β) and φ 2 (α, β) independent of both γ and δ defined as Equation (5) refers to the measure of fat-tail behaviors with the focus on estimating α, and equation (6) is a measure of skewness effects with the focus on estimating β. With empirical values of sample quantiles and employing linear interpolation with tabular lookups, the estimatesα,β are inversely obtained. To avoid α being larger than 2, outside the parameter range, can be no larger than the upper range 2.439, which corresponds to the case of α = 2 (note that β is not identified in this case).
Next, the scale and location parameter γ and δ can be estimated using the functions defined as The function φ 3 (α, β) indicates the standardized form of sample sizes for the middle part of distribution. Since it does not depend on γ nor δ, the value can be informed by tabular look-ups based on α and β, which the relations are studied and known beforehand. After calculatingγ =x 0.75−x0.25 in equation (7), the location parameter δ can be estimated from equation (8) using the values φ 4 (α,β) andγ. The relations of the parameter values and the function value φ 4 (α, β) are again, studied and known beforehand. In the case of α = 1, φ 4 (α, β) diverges and we cannot obtain the estimates for δ. McCulloch therefore suggests a complicated approach to overcome the discontinuity of the stable CF. The method improves other issues and provides accurate estimates, however, it has parameter restrictions and can be applied only when α ≥ 0.6.

B. Characteristic function-based method
The CF-based method relies on the use of a consistent estimator of the CF ϕ(k) for any fixed k. The advantage of this method essentially lies in the fact that the stable CF can be expressed explicitly, making discussions straightforward compared to methods based on other distribution forms. Under the assumption that given data X n (n = 1, 2, . . . , N ) are ergodic [23], the CF is obtained empirically by the following equation, There are several approaches for estimating parameters of stable laws that take advantage of the explicit form of CF. Koutrouvelis (1980) [14] proposed a regression-type approach, which employs the iteration of two regression runs. Moreover, the regression of the method requires different values of initial points k depending on initial estimates of the parameters and sample sizes. The number of points necessary for the regression also varies over initial conditions. Although the accuracy of β is unsatisfactory in some cases, the method generally shows accurate estimates of α, and hence it is often suggested as a practical method for empirical analysis [24,25]. However, some studies compare the method to McCulloch's quantile method and report that the regressiontype method does not significantly improve the classical quantile method [26,27], especially for α smaller than 1.
Other studies simplified the method by eliminating the iteration process and fixing the initial points to some extent, but still leaves behind the issues of estimating when α is small [28,29]. We do not consider the regression-type approach in this paper as the method generally relies on iteration and the estimates cannot be written analytically.
Another approach is based on the method of moment [15], which was later remodeled and simplified with the use of two given points of the CF [16][17][18]. Starting off with the CF with the points k 0 and k 1 , taking the absolute value cancels out the effect of parameters β and δ, and we obtain Taking the cumulant function, which is the natural logarithm of the CF, leads to the same discussion neutralizing the effect of parameters β and δ. The equation ln ϕ = ln |ϕ| + j(arg ϕ + 2nπ) implies that the real part of the cumulant function corresponds to the natural logarithm of the absolute value of CF, shown as for any value of k. We consider only the positive values for convenience, since the CF is a symmetric function. By solving the above equations simultaneously, parameters α and γ can be estimated shown aŝ Since the one-parameterization form in equation (2) is discontinuous at α = 1, the estimation of the remaining parameters β and δ is divided into two cases. When α = 1, the cumulant function of stable distributions with the points k 0 , k 1 > 0 are As we need the information of the parameters β and δ, we take the imaginary part. Then the parameters β and δ are estimated by solving the above equations simultaneously and using the estimatesα andγ: In the case of α = 1, the CF takes a discontinuous form and the cumulant functions are written as Then the parameters are estimated by solving the above equations simultaneously as well: For simplicity, we express the estimates as a function of given points k 0 and k 1 as follows: whereβ andδ additionally needs the information of the estimatesα andγ. Sometimes, the estimates can possibly outrange the parameter spaces α ∈ (0, 2], β ∈ [−1, 1], and γ > 0, especially when the true parameters are close to the borders. In such cases, the parameters are set to the closet border, except for α and γ, the estimates are set no lower than 0.01. Applications with other parameterizations use slightly different forms of CF, but the stable parameters are estimated essentially by the same procedure as explained above. For the zero-parameterization, which is another common parameterization form, the CF is replaced to its corresponding form shown in equations (3) and (4) for equations (10) and (14) (or (17) Next α is estimated as shown in equation (12). Then, the scale parameter in our criterion,γ, is obtained as, The remaining parameters β and δ are then estimated straightforwardly as similar to the case of the oneparameterization form. By replacingγα toĉ in equations (15) and (16) (or equations (18) and (19)), and using the points k 0 and k 1 give the estimates.

IV. PROPOSED APPROACH FOR THE CHARACTERISTIC FUNCTION-BASED METHOD
In this section, we make an improvement of the CFbased method by discussing how the points related to the CF should be chosen. We propose a technique that provides a flexible selection of the points. We also clarify the difference of how the points are selected between our proposal and the procedures in other existing CF-based methods.
A. Argument for the inference of point k1 Two positive points of the CF, k 0 and k 1 (k 0 = k 1 ), are ought to be selected to identify all four parameter estimates. As mentioned before in this paper, the absolute value of the CF in equation (10) is independent of the skew and location parameters for any k, and provides information of α and γ. When k = 1/γ is satisfied, the absolute value of the CF takes a constant value The advantage of setting k = 1/γ as one of the candidate points is to reduce any estimation bias influenced by certain parameter values since we expect to get a constant estimate which is independent of all four parameters. When γ 1, however, empirically obtained values can cause significant estimation errors for the scale parameter in equation (26) [18,32]. Therefore, we first consider a temporary estimate of the scaling parameter, γ, just in case the data exhibits scale far from the standardized form (γ = 1).
Take the natural logarithm of equation (26). The temporary estimate can be obtained by approximately solving the equation that numerically satisfies using a simple one-dimensional search function [33], or any other optimization procedure. Our rough estimatẽ γ is then used for standardizing, or pre-standardizing, the candidate points. Specifically, point k 1 is set to 1/γ, where ln |ϕ(k 1 )| empirically takes −1.
As explained above, pre-standardization is preferred especially when we suspect that datasets have too large or small scales. Whenever a new set of points is required for the parameter estimation process, we conduct prestandardization. Point k 1 is replaced to 1/γ, whereγ is the latest scaling parameter estimate available at that time.

B. Argument for the inference of point k0
For the argument of selecting point k 0 > 0, which is perhaps the most important proposal in our study. We focus on the absolute value of the CF. Bibalan et al. (2017) proposed to calculate the distance between two absolute values of CFs with different index parameters α, the Gaussian case (α = 2) and the Cauchy case (α = 1) [16]. They set k 0 > 0 to the point which corresponds to the maximum distance and the other point to k 1 = 1. Although the absolute CF changes depending on the index parameter α, their approach considers a fixed distance and essentially chooses an identical point for any value of α ∈ (1, 2]. In addition, the distance they consider does not account for the case of α ∈ (0, 1]. Our approach is an extension of Bibalan et al. (2017), and provides a more generalized technique of selecting the points. We deal with the problem that the distance between two absolute values of CFs can vary depending on the parameters. The basic idea is to find the point where the absolute CF, |ϕ(k; α)|, presents the maximum sensitivity with respect to α. In other words, we discuss the point where the distance between the absolute CF of index parameter α, |ϕ(k; α, β, γ, δ)|, and the absolute CF of α + ∆α, |ϕ(k; α + ∆α, β, γ, δ)|, shows the largest distance. Such a point is considered as k 0 in our study.
To make our discussion more simple, we consider the absolute CF as a function of variable η: where η = γk (k > 0) is a newly introduced variable which depends on γ and k. The distance can be expressed as exp(−η α+∆α ) − exp(−η α ) . The candidate point for η 0 = γk 0 , where the maximum distance is achieved, can be calculated by Solving this equation for η > 0 yields two solutions, η ∈ (0, 1/γ) and η ∈ (1/γ, ∞). For both points, the absolute value of CF shows the largest ratio of change in a local sense. The smaller point η ∈ (0, 1/γ) is employed, because the distance at the smaller point tends to have larger values than that at the larger point η ∈ (1/γ, ∞), which enables us to estimate α and γ in a more desirable and informative manner. Another reason is that smaller |k| is preferred rather than larger |k|. As |k| → 0, the asymptotic variance of the empirical cumulant function decreases [18]. With empirical CF obtained by i.i.d. distributed datasets, the relation holds [34], which implies that as k becomes larger, the empirical absolute CF |ϕ N (k)| is likely to be subject to sample errors. Thus, the smaller η = γk should be considered in this study.
The above discussion implies that k should be set close to zero (but not at zero because then the CF takes a constant value and no information of the parameters will be provided). But at the same time, the employed smaller point is standapart from zero to some extent, so that the empirical CF will be more or less exposed by sample errors. Therefore, the choice of points derived from equation (28) is unsatisfactory, and hence the distance exp(−η α+∆α ) − exp(−η α ) should be modified. To reduce the effect of sample errors, we introduce a weight function w(η) that decreases monotonically as η becomes larger (note that the introduced variable η = γk has a linear relationship with k).
This choice leads to the association of the weighted distance with a statistical measure used for goodness-of-fit tests, developed by Matsui and Takemura (2008) [35]. They propose the following test statistic based on empirical CFs, where h(t) is a monotonically decreasing weight function. D N,κ denotes the weighted L 2 -distance between the empirical CF and the symmetric standardized stable CF ϕ(t; α, 0, 1, 0). This weighted L 2 -distance can be associated with the weighted distance we are considering now. Taking the absolute value of a CF yields again a standardized form of a CF with β = 0 and δ = 0: exp(−η α ) = |ϕ(k; α, β, γ, δ)| = ϕ(η; α, 0, 1, 0). Thus, the absolute values of CF with index parameter α and α+∆α are equivalent to the symmetric standardized stable CFs, ϕ(η; α, 0, 1, 0) and ϕ(η; α + ∆α, 0, 1, 0), respectively. The weighted L 2 -distance between these CFs essentially coincides D N,κ , when the weight function satisfies for η > 0. In this case, the difference between the CFs can be evaluated more accurately with the background of a meaningful measurement. Following Matsui and Takemura (2008), the asymptotic distribution of D N,κ is numerically evaluated and the critical values of the test statistics are approximately obtained [35]. Through computational simulation, they provide evidence that the test is most powerful when κ = 5.0 (h(η) = exp(−5|η|)), especially for heavy tailed distributions. Thus, our choice of the weight function is w(η) = exp(−2.5|η|), since τ = κ/2. Other weight functions such as w(η) = exp (−|η|) and w(η) = exp (−η 2 ) [32,36] can be employed, but lacks a conclusive evidence for the use of these alternatives.
With the weight function, the candidate points η > 0 are calculated by solving the following equation: where τ = 2.5. Then we have  (31), is shown in the solid black line. The blue plot shows the simulated results for the best point with the minimum MSE for α and β over 100 simulations. We consider the MSE of α + 1 10 β because the accuracy of β is generally worse roughly by ten times than the accuracy of α, and also that β estimates are usually susceptible to α estimates [9,14,16]. The simulation is implemented for each value of α ranging within the parameter space of 0.2 to 1.95.
the maximum rate of a change, or the maximum sensitivity, of the absolute CF with respect to α.
There could exist some relationship between α and η since they are interrelated due to g(α, η) = 0. When some estimateα is given, the corresponding point is obtained by computing η that satisfies g(α, η) = 0, and vice versa (the corresponding parameter α of a given pointη can be calculated by computing the equation g(α,η) = 0). As we have discussed previously in this subsection, we focus on the point closer (smaller) to zero out of the two candidates of the calculated points from equation (31). Figure 2 ascertains whether our approach of equation (31) correctly estimates the parameters of stable distribution. The model clearly characterizes the distinctive relationship between α and η, which are empirically verified via simulation using synthetic data generated from random stable variables [37]. This indicates that our selection of points is valid for identifying desired points in the estimation process.
In practice, α is unknown. Hence the selection of point η 0 = γk 0 is undecidable, so that the parameters for the stable law cannot be estimated directly. To cope with this problem, we first aim to get a rough estimate of α calculated by using the temporary scale estimateγ. The rough estimate is considered poor as the estimation method, but it plays a role in starting off the estimation process with reasonable initial values. The accuracy of both points (η 0 = γk 0 and η 1 = γk 1 ) and the parameters (α, β, γ, δ) can be improved by alternating searches of α and η from our relation model g(α, γ) = 0 several times to get sophisticated estimates. With estimates η 0 and η 1 , the four parameters are ultimately calculated.

C. Estimation procedures
Here we present our proposed algorithm for the estimation of all four parameters of stable laws by utilizing the relationship between α and η. Regarding the fact that empirically obtained estimates occur substantial errors induced by γ 1, we conduct a pre-standardization with k replaced to η = γk. Using the expressions of the estimates in equations (20) (21) (22) (23), our algorithm is written as follows: Step 1: Compute a temporary estimateγ temp from sample data X n (n = 1, 2, . . . , N ) that satisfies the equation, Step 2: Set where ξ is any initial value of ξ ∈ (0, 1).

V. NUMERICAL ASSESSMENTS
In this section, we show numerical assessments for the estimation of stable laws. We compare the performances of our proposal approach to other state-of-art approaches using the MSE and the KS-distance. The comparison is studied for three approaches. We focus on the approaches of characteristic function-based methods presented by Bibalan et al. (2017) [16] and Krutto (2018) [18]. We also compare with the traditional QM method [8,9] explained in subsection 3.1, to provide a benchmark with a well-known criterion. Note that all three approaches above exhibit closed-form expressions for all four estimates of stable parameters. Bibalan et al. (2017) have shown that their approach generally outperforms other methods that yield a closedform expression, such as the FLOM, the QM, and the MOLC [16]. Krutto (2018) also compares the performances with several well-known methods and concludes that the method gives accurate estimates [18]. Since both of them belong to the family of the CF-based method, the selection of the points k 0 and k 1 plays an important role. In Bibalan et al. (2017), k 1 is set to 1. Point k 0 is always set to where the point shows the maximum distance between the absolute Gaussian CF and the absolute Cauchy CF, by using the estimates of γ α which they are calculated beforehand. It should be mentioned that the CF in this case poses a alternative definition of the scaling parameter, so we eventually obtain γ in the last procedure in equation (25). On the other hand, Krutto (2018) suggests to employ two points that satisfies ln |φ(k 0 )| = −0.1 ln |φ(k 1 )| = −0.5, under empirical searches [18]. We examine the performance for each parameter of stable distribution in addition to the fit with the entire estimated stable distribution. We also refer to the effects of sample sizes for each estimation method. For all the simulations in this paper, we generate L = 500 synthetic data of N = 10000 i.i.d. random stable samples. Synthetic random data sequences following a stable distribution can be generated by algorithms constructed by Chambers et al. (1976) [38], Weron (1996) [37], and Umeno (1998) [39]. Umeno (1998) generates random stable variables based on the superposition of chaotic processes. The classical method of Chambers et al. (1976) is widely known as the pioneer of all the methods, which the algorithm was reorganized and corrected later by Weron (1996). Weron's algorithm is our choice of method, which is simple and is the fastest in calculation.

A. Performance of parameter estimates
The performance of the estimated parameters are examined by the MSE criterion: where θ and L = 500 is the parameter of stable laws and the number of times the simulation is implemented, respectively. We calculate the MSE of all four parameters and evaluate each parameter individually. Table I shows the simulation results of the MSE associated with the estimate bias for each parameter. We consider the cases of parameters with α = 0.5, 1.5, 1.8 and β = 0, 0.5, all with a standardized form of γ = 1 and δ = 0. Our proposed approach generally provides the most accurate estimation with the smallest MSE. Especially for the index parameter α and δ, our approach significantly improves the accuracy of the estimates. Note that for the QM, the method has parameter restrictions of α ≥ 0.6 and hence the cases with α smaller than 0.6 can not be implemented. More detailed simulation results for different cases of parameters are shown in Figures 5, 6, 7 in Appendix A. In particular, we show the cases of S(α, −0.1, 1, 0), S(0.8, β, 1, 0) and S(1.6, β, 1, 0), with parameter values varying within the parameter ranges. The results imply that for whatever parameter combination, our method generally outperforms the others with the highest accuracy. Although we find that other methods sometimes show higher accuracy on either the parameter α or δ in cases of 0.6 ≤ α ≤ 1.2 in S(α, −0.1, 1, 0) in Figure 5 and −0.3 ≤ β ≤ 0.3 for S(0.8, β, 1, 0) in Figure 6, our method appears to be powerful for estimating all four stable parameters.

B. Performance of the estimated distribution
Next, we examine the performance of estimating stable laws from a different perspective; evaluation of the entire distribution. We use the KS distance expressed as which represents the maximum distance between two distributions in terms of the CDF. Here P (x) andP (x) denotes the empirically obtained CDF, and the theoretical estimated CDF, respectively. The standard density and distribution functions of stable distributions are numerically derived approximately by implementing the Fourier integral formulas [40,41], which are available in package libstable that provides good approximation values [42]. KS distance is one of the most major standards for numerical assessments when discussing stable laws. We set aside any issues related to numerical approximations of stable distributions, so that we can focus on the performance between the methods. The root mean square (RMS) of the KS distance is used for the numerical assessment to make the small differences of the comparison results more apparent.  1, δ). The RMS of the KS distance is calculated for each case with various values of parameters ranging within parameter ranges. We find in Figure 3 (c) that the estimation for the scaling parameter γ = 1 poses significant estimation errors. This is caused by the effect of sample errors induced by the scaling parameter γ far from the standardized form, as shown in equation (26).
On the other hand, our proposed method achieves the smallest value of KS distances for all cases of parameter combinations. This proves that we are also successful in improving the estimation of the entire stable distribution.     Needless to say, the accuracy of the estimation method strongly depends on the number of samples. Larger sample sizes give more information of the dataset whereas smaller sample sizes have only little information making it challenging to detect the true values. We examine the effect of sample size by comparing the performance among the estimation methods. Figure 4 displays the MSE of each parameter of stable distribution as the sample size N changes from 300 to 10000. The study is examined for the case of S(1.4, 0.2, 1, 0). The MSE simulated by means of our method decreases with the order O(1/N ) while the MSE simulated by means of other representative methods also exhibited similar behaviors of order. Our proposed approach offers the best performance except for the location parameter δ, where the QM method sometimes give more accurate estimates for large datasets.

VI. APPLICATION TO FINANCIAL EMPIRICAL DATA
This section shows application of the proposed estimation method to real financial data. We provide several empirical studies to present that our proposed approach is appliable for a wide range of empirical analysis in finance.
Asset price returns in various financial markets tend to show interesting properties of stable laws ever since Mandelbrot (1963) first revealed that stable distribution fits cotton price returns better than the classical Gaussian distribution [43]. This argument have attracted attention to identifying price behaviors in many financial fields such as equities [44][45][46], price consumer index inflation [47], metal markets [48], oil markets [49], and Cryp- tocurrency markets [34]. We investigate return distributions of the Japanese Yen currency exchange rate in terms of the US dollar (USDJPY) and the West Texas Intermediate (WTI) crude oil futures market, both of which are potent indices in finance. The basic statistics of the indices are provided in Table II. We explore both cases of common daily analysis and high-frequency data analysis. In particular, we use daily and one-hour return time series for the USDJPY and the WTI market, respectively. Since the scale of returns for both cases are too small for the method based on Bibalan et al. (2017) to give plausible estimates, we do a pre-standardization process beforehand. We multiply returns by 100 and after the estimation the parameters γ and δ are adjusted by dividing them by 100. Table III presents the estimates of the fitted stable distribution associated with the KS-distance between the empirical distribution and the estimated stable distribution for USDJPY, calculated based on four controversial estimation methods. Our primary focus is on the KS-distance value. The results show that the estimated distribution based on our proposed method presents the smallest value among other estimation methods. The smallest KS-distance implies that our method exhibits stable laws that best describes the observed data. Parameter estimates and the distance measure for the WTI market are shown in Table IV. The result indicates that the outstanding performance of our method also holds for high-frequency data with the lowest KS-distance. What makes the development of the estimation method a crucial matter is that the parameter estimates can differ so much among the methods when applied to empirically observed data, even for large datasets. We find in Table IV

VII. CONCLUSION
This paper has proposed a new approach for estimating stable laws and applied this approach to the exploration of price behaviors in financial markets. Our new technique is developed under the method of moments, which is one of the widely known CF-based methods that require the choice of appropriate momental points. The points necessary for the estimation process are flexibly chosen, as the estimation accuracy of stable laws depends heavily on their true parameter values. We have focused on the fact that the index parameter α and the desired momental points exhibit a distinctive relationship, which is a new perspective in the literature. This relation is modelled as g(α, η) = 0, based on the idea of employing points η at which the weighted absolute values of the CF present the maximum sensitivity. To detect appropriate points, we have suggested a procedure relying on the combination of empirical searches and algorithmic approaches. The advantage of employing these points is that the parameters of stable laws can be estimated in a more precise manner while remaining straightforwardly the implementation of the method. The relative performance of the parameter estimates is benchmarked against other existing methods, specifically the QM and the methods of Bibalan et al. (2017) and Krutto (2018), through simulation studies in terms of the MSE and KS-distance criteria. The results have implied that our method is the most powerful with the best performance. Our approach assures that the estimates of all four parameters of stable laws present a closed-form expression without any restrictions on parameter ranges, making the method significantly practical. We have also explored the behaviors of price fluctuations in several financial markets to show that our method is applicable for empirical financial studies. For the USD-JPY exchange rate and the WTI crude oil future price, our method supports stable laws with the highest performance among all the other methods discussed in this paper. This would motivate us to further develop analytical methods for examining stable laws, as well as to further investigate various features of financial markets.  Comparison of the MSEs for the methods based on the proposed approach, Bibalan et al.'s approach, Krutto's approach, and the QM method. The MSEs of each stable parameter are studied for cases of parameters α = 0.8, γ = 1, and δ = 0 with β ranging from -0.9 to 0.9 (N = 10000, L = 500). Comparison of the MSEs for the methods based on the proposed approach, Bibalan et al.'s approach, Krutto's approach, and the QM method. The MSEs of each stable parameter are studied for cases of parameters α = 1.6, γ = 1, and δ = 0 with β ranging from -0.9 to 0.9 (N = 10000, L = 500).