{'fraction_non_alphanumeric': 0.0447255369928401, 'fraction_numerical': 0.012887828162291169, 'mean_word_length': 4.475953998954521, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Analysis and application of multiplicative stochastic process with a sample-dependent lower bound

A multiplicative stochastic process with the lower bound lognormally distributed is investigated. For the process, the model is constructed, and its distribution function (involving four parameters) and the related statistical properties are derived. By adjusting the parameters, it is confirmed that the theoretical distribution is consistent with empirical distributions of some real data.


I. INTRODUCTION
Probability distributions provide useful information on several complex phenomena and help achieve effective modeling. Even if the precise dynamics or underlying mechanisms are unknown, a probability distribution can be simply constructed by observed data. As such, the analysis of the empirical probability distributions is regarded as a possible first step toward understanding complex systems.
The heavy-tailed probability distribution, which decays slower than exponential distributions, is an important class of probability distributions [1]. Power-law distribution [2] is representative of the heavy-tailed distributions. The appearance of a heavy-tailed distribution implies that there exist elements having a far larger size than the characteristic size and that a mechanism for generating large-sized elements is embedded in the system.
Meanwhile, lognormal distribution is another representative heavy-tailed distribution that appears in various complex systems [3,4]. Here, the random variable X is said to follow a lognormal distribution if its logarithm ln X is normally distributed. The probability density f (x) and cumulative distribution F (x) of the lognormal distribution can be described as follows: where µ and σ correspond to the mean and standard deviation of ln X, respectively, and erf is the Gauss error function defined by the following: erf(x) = 2 √ π x 0 e −y 2 dy.
The lognormal distribution is generated by the following multiplicative stochastic process where B 0 , B 1 , . . . represent positive random variables that are independently and identically distributed. In fact, according to the central limit theorem, the distribution of ln X n converges to a normal distribution. Qualitatively, because the lognormal distribution is attained via numerous factors that are accumulated multiplicatively, the lognormal distribution is observed in a variety of situations. This model was first introduced by Kolmogorov [5] to explain a fragment-size distribution in a cascade fracture. The multiplicative stochastic process (2) can produce a power-law distribution when an additional effect is introduced, one example being the introduction of a lower bound a for X n [6] such that X n is restricted to not less than a: In the n → ∞ limit, X n follows a stationary power-law distribution: where the exponent β > 0 is given by the positive solution of The corresponding probability density function is expressed as Historically, this model dates back to the work of Champernowne [7], who applied it to the power-law income distribution. A continuous version that involves geometric Brownian motion with a reflection boundary has been investigated [8]. Other mechanisms for the power-law distributions from the multiplicative process (2) have been also studied in terms of the introduction of additive noise [9], a reset event [10], random stopping [11], and the sum of variables [12].
In this paper, we investigate a stochastic process based on business-related empirical data. Here, the empirical distribution of the data follows a lognormal distribution up to the mid-sized region, with a tail fatter than the usual lognormal distribution. We propose that this tail behavior can be modeled via a multiplicative stochastic process with a specific lower bound, in which the initial value is drawn from a lognormal distribution and the lower bound is equal to this initial value. After formulating the model, we derive the distribution function and its properties.
II. MODELING Figure 1 shows the cumulative number of troubleshooting costs in a specific Japanese company. The cumulative number is defined as the number of cases with costs larger than x. The absolute values are substituted with the costs divided by the average as shown in Fig. 1 [13]. The total number of cases is N = 2571. The solid curve corresponds to the lognormal distribution with µ = −0.84 and σ = 1.24. (Note that the cumulative number is given by N F (x).) The empirical distribution deviates from the lognormal distribution in the tail (x 10).
While only a few (the order of 10 1 ) points deviate from the lognormal distribution, they are highly expensive cases, ranging from approximately 10 to 100 times greater than the average value. Hence, the precise description of the tail behavior based on a mathematical mechanism will be financially important in practical situations. Moreover, it has been widely observed that the main part of the data follows a lognormal distribution, but the tail deviates from the lognormal [14][15][16]. Such tail behavior is generally difficult to determine and can often be controversial. For example, while the population distribution of American cities was reported to be lognormally distributed [17], its tail was later noted to exhibit a power-law decay [18]. This issue has not been resolved yet because the measurement of the tail is problematic [19,20]. Furthermore, the population distribution varies according to city, town, and village [21], and a numerical model incorporating the natural and social effects on population increase or decrease in a simplified form has been proposed [21,22]. Therefore, the tail of an empirical distribution is worth modeling and analyzing in theoretical terms.
Illustration of the interchange of the integration in Eq. (7). The shaded area composed of semi-infinite triangular (corresponding to G1(x)) and rectangular (G2(x)) areas represent the domain of the iterated integral in the y-a plane.
Assuming the tail of the empirical distribution shown in Fig. 1 follows power-law decay, we obtain the power-law exponent −1.53 ± 0.08 using the method by Clauset, Shalizi, and Newman [23]. However, this power-law decay (shown by the dashed line in Fig. 1) does not seem to express the tail of the empirical distribution well. The power-law decay cannot always be justified theoretically, and it should be regarded as a rough assumption or a working hypothesis. A more appropriate function may be required to understand the system accurately.
We thus propose the following phenomenological model. Based on the lognormal distribution in the low-and midcost regions, it is assumed that the cost is determined by several multiplicatively affecting factors. Thus, this can be naively modeled using the stochastic process described in Eq. (2), which produces the lognormal distribution. Additionally, for expensive cases, additional tasks, such as individual inspection and customization, tend to be involved. These additional tasks may increase the cost but never reduce it to less than the initial cost. This operation can be expressed by a multiplicative process wherein the initial cost acts as the lower bound. Due to this additional process, the price can further increase, resulting in a deviation from the lognormal distribution.
The above model can be simulated as follows. First, let us define x * as the boundary cost where the empirical data begins to deviate from the lognormal distribution. A random number x is then drawn from the lognormal distribution. If x ≥ x * , the following multiplicative process, wherein the lower bound is x itself, is performed, where the initial value is X 0 = x and the evolution is given by until X n becomes stationary. Here, B 0 , B 1 , . . . are independently and identically distributed positive random numbers. Unlike in Eq. (3), where the lower bound a is constant, the lower bound x is distributed lognormally and varies for each sample in this model. Therefore, x is lognormally distributed for x < x * , whereas the distribution changes for x ≥ x * due to the additional process. In the next section, we analyze the proposed model and derive the distribution function for x ≥ x * .

III. MODEL ANALYSIS
As mentioned earlier, the proposed model involves a multiplicative process with a varying lower bound. First, we derive the corresponding probability density g(x) and cumulative distribution G(x).
For x < x * , the effect of the lower bound does not come into play. Hence, g(x) and G(x) are identical to the lognormal distribution f (x) and F (x) described in Eq. (1), respectively. For is obtained by superposing the stationary distributions corresponding to the different values of the lower bound. When the lower bound is a, the stationary distribution follows the power law φ a (x) described in Eq. (5). In the proposed model, the lower bound a has a lognormal weight f (a) [see Eq. (1a)], such that The cumulative distribution can be expressed as follows: The order of integration can be interchanged as follows (see Fig. 2 for reference): By the normalization ∞ a φ a (y)dy = 1, we can easily obtain The function G 2 (x) can be written as Since the integrand represents the lognormal probability density (as a function of a) where the parameter µ is replaced with µ + βσ 2 , we have Thus, for x ≥ x * , we can obtain Here, the first term on the right-hand side represents the lognormal distribution, while the second term is regarded as the attendant correction term. By differentiating Eq. (11), we obtain Figure 3 shows graphs of g(x) and G(x) with fixed (µ, σ) = (0, 1.5) and x * = 10, and varying β = 0.5, 1, 2, 4, and 8. The tail decays rapidly for large β, and the probability density g(x) is not continuous at x = x * : g(x) → 0 as x ց x * . In accordance with the discontinuity of g(x), the G(x) graph is not smooth at x = x * .
If a random variable X has the probability density g(x), the kth moment of X can be calculated as  x * = 4.21 Ratio Normalized cost x for k < β, and the moment diverges for k ≥ β. The calculation of the integral is straightforward but cumbersome. This criterion for the divergence of moments can be derived simply using the asymptotic power law g(

IV. DATA ANALYSIS
The function G(x) involves four adjustable parameters: µ, σ, x * , and β. The optimal value of these parameters for fitting the given data can be determined as follows. For fixed x * , we can obtain µ and σ according to the lognormal fitting in the region x < x * . Using these values µ and σ, we can obtain β by fitting Eq. (11) in the region x > x * . By calculating µ, σ, β, and the fitting error for different x * , we can determine the optimal set of the parameters as one minimizing the fitting error. Figure 4 shows the result of the curve fitting of G(x) to the empirical distribution shown in Fig. 1. Here, the parameter values are (µ, σ) = (−0.84, 1.24) (same as those in Fig. 1), x * = 4.21, and β = 5.28. G(x) is consistent with the empirical distribution almost over the entire range, which indicates the validity of the proposed model i.e., the multiplicative stochastic process with the sample-dependent lower bound.
The naive power-law fitting shown in Fig. 1 gives the power-law exponent −1.53(> −2) for a cumulative distribution, implying that the variance diverges if the same power law extends to infinity. In contrast, we obtained β = 5.28 when using G(x), which results in finite variance. Here, it can be expected to obtain the correct statistical properties of the data when using G(x), ensuring appropriate and effective cost management in practical situations.
As outlined in Sec. III [see Eqs. (8) and (10)], G(x) for x ≥ x * can be decomposed into two parts G(x) = G 1 (x) + G 2 (x), where G 1 (x) is a lognormal function, and G 2 (x) is a power-law-like function (not an exact power law due to the erf factor). Figure 5(a) shows G 1 (x) and G 2 (x) for (µ, σ) = (−0.84, 1.24), x * = 4.21, and β = 5.28 (same values as in Fig. 4). In this case, the magnitude of G 2 (x) is greater than that of G 1 (x) for x > 25.8, and the power-law-like part G 2 is dominant in the tail. Specifically, as Fig. 5(b) shows, the ratio G 2 (x)/G(x) represents an increasing function. Furthermore, using the asymptotic expansion 1 − erf(z) ≃ exp(−z 2 )/(z √ π), G 2 (x)/G(x) → 1 can be presented as x → ∞ for any parameter value. In addition to the troubleshooting costs, we present the empirical distributions of other social data. Figure 6(a) shows the distribution of the number of employees aggregated according to a municipality (N = 1737) in Japan [24] in the form of the cumulative number. Here, the solid curve represents G(x) for (µ, σ) = (9.17, 1.59), x * = 3.09 × 10 5 , and β = 4.20. We believe that the lognormal behavior in the x < x * region is related to the lognormal population distribution [21], and that the heavy tail likely reflects the extreme concentration of companies in metropolises, such as Tokyo and Osaka. Meanwhile, Fig. 6(b) shows the cumulative number (points) of the management expenses grants for 90 national universities in Japan in 2019 [25] and G(x) for (µ, σ) = (2.11, 0.82), x * = 21 × 10 9 yen, and β = 3.80 (solid curve). Here, 12 universities have x greater than x * . While it is not clear whether these data are generated by the same statistical mechanism as that outlined in this study in terms of the troubleshooting costs, we believe that G(x) effectively expresses the empirical distribution; thus, providing some insight into these phenomena.
Next, we compare the proposed distribution G(x) with other distribution functions. We selected the following three distributions that have been frequently used in the analysis of empirical data: the lognormal distribution F (x) [shown in Eq. (1b)], the lognormal distribution with power-law tail, and the double lognormal distribution (a mixture of two lognormal distributions) given by where µ 1 , σ 1 , µ 2 , σ 2 , and α are the fitting parameters. Figure 7 shows the cumulative number of employees in municipalities (a) and management expenses grants for the Japanese national universities (b) (the datasets are the same as with Fig. 6), and these three distributions.
To obtain the center graphs of Fig. 7, we performed the following procedure. First, the power-law exponent and the lower bound of the power-law distribution were determined using the method outlined in Ref. [23], and we regarded this lower bound as the crossover point of the power-law and lognormal distributions. Next, the lognormal fitting was performed in the range up to the crossover point. Finally, the coefficient of the power-law distribution was computed so that the power-law and lognormal functions were continuously connected.
To evaluate the quality of curve fitting using the four distributions (the lognormal distribution, the lognormal distribution with power-law tail, the double lognormal distribution, and the proposed distribution), we measured the distance between the empirical distribution and each candidate distribution. Among a wide variety of distances or metrics between probability distributions, we selected two simple ones: the Kolmogorov-Smirnov statistic D and the first Wasserstein distance W 1 . If the cumulative distributions of the empirical data and the fitting result are respectively written as F emp and F fit , these two distances are defined as These distances take non-negative values, and the small distances indicate that the estimated distribution is close to the empirical one.
The calculation result of W 1 and D for the lognormal distribution, lognormal distribution with power-law tail, and double lognormal distribution in Fig. 7, in comparison with the proposed distribution G in Figs. 6 and 4, is presented in   Table I. The relative values (the distance between each distribution and the empirical distribution divided by the case of the proposed distribution G) are shown in the table; the value greater than unity represents that the distribution G is closer to the empirical distribution than the corresponding distribution with respect to the corresponding distance.
In the data analysis of complex systems including social ones, the tail of empirical distributions has often been focused on. However, the above statistical distances are not suitable for the evaluation of the tail behavior. In fact, Fig. 8 Fig. 7(a)]. The differences in the tail (x 10 5 ) are small compared with those in the body (x 10 5 ). and similar W 1 and D values in Table I are considered to be attributed to this common lognormality. To magnify the difference in the tail, we propose two distances: where x min and x max are the minimum and maximum values of the empirical data, respectively. (The upper limit of the integral is set to x max in order to avoid the singularity F emp (x) = 0 for x > x max .) Like other statistical distances, ∆ 1 and ∆ 2 are always nonnegative and they become smaller as the estimated distribution F fit is closer to the empirical one F emp . The quality of curve fitting based on ∆ 1 and ∆ 2 is presented in Table II. As in Table I, the relative values of the three distributions in Fig. 7 against the distribution G in Figs. 6 and 4 are shown in this table. The proposed distribution G describes the tail of empirical distributions better than (or at least comparable to) the other three distributions (lognormal, lognormal with power-law tail, and double lognormal).
In addition to the quantitative aspect, another advantage of G is that it possesses a mechanism (stochastic process) different from the other three distributions. That is, the agreement of the empirical dataset with the proposed distribution implies that the data was generated under the effect of the proposed mechanism, namely a multiplicative stochastic process with a sample-dependent lower bound. This conclusion is based on only the three specific datasets shown in Figs. 4 and 6, and general discussion on the quality of each candidate distribution is the subject for a future study.

V. DISCUSSION
The empirical data of troubleshooting costs outlined in this paper is normalized in terms of the average. Here, we analyze the effect of this normalization. If we set m as the average cost, the raw valuex and the normalized value x satisfy x =x/m. The threshold x * is transformed intox * = mx * , since x ≶ x * changes tox ≶ mx * =x * , while the parameter µ in G(x) appears in the form of ln x − µ or ln x * − µ at all times [the factor exp(µβ)x −β is reduced to exp(−β(ln x − µ))]. Therefore, the parameter µ becomesμ = µ + ln m. It should be noted that the parameters σ and β in G do not change due to the normalization. In conclusion, writing the parameters explicitly as G(x; µ, σ, x * , β), we have G(x; µ, σ, x * , β) = G(x; µ + ln m, σ, mx * , β). The same dependence of the parameters holds even if the constant m is not equal to the average value.
The distributions g(x) and G(x) are derived based on the assumption that the lower bound of the multiplicative process is lognormally distributed. We can generalize these results by replacing the lognormal weight f (a) described in Eq. (6) with alternative probability distributions. Here, when a probability density h and its cumulative distribution H are used instead of the lognormal distribution f and F , we obtain G 1 (x) = H(x) and which are similar to Eqs. (8) and (9). The asymptotic power law G(x) ≃ G 2 (x) ∼ x −β generally holds if h decays faster than x −β−1 .
The tail of an empirical probability distribution usually contains only a small amount of data. For example, while the total number of cases is N = 2571 in Fig. 4, only 85 cases (3.3 % of the total) lie in x ≥ x * . Due to statistical fluctuations, the distribution tail is difficult to specify. From Figs. 4 and 6, the proposed distribution G fits well even up to the tail, which appears to increase reliability for the estimation of the asymptotic power-law exponent β. By its construction, the proposed model, a multiplicative process with a sample-dependent lower bound, will be useful for a lognormal distribution with a fat tail, especially for a system whose dynamics or evolution rules change at a certain threshold value. In fact, the troubleshooting cost, the number of employees, and the university management expenses grants analyzed in this study are real examples which follow the proposed distribution G, and this result suggests that the proposed mechanism may be embedded in the generation process of them. Moreover, we expect that there exist other systems that can be well described by the proposed model.
We think that the proposed distribution g or G is an effective distribution but not the definitive one. In general, the distribution function of a given dataset should be determined by examining several candidate distributions, and we hope that the distribution G outlined in this paper becomes one of the possible choices.