Harmonic distribution

Harmonic

Probability density function ProbDensFunc
Cumulative distribution function CumDisFunc
Notation	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mathrm{Harm}(m,a)\,
Parameters	m ≥ 0, a ≥ 0
Support	x > 0
PDF	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \frac{1}{2xK_{0}(a)}\exp\left(-\frac{a}{2} \left(\frac{x}{m}+\frac{m}{x} \right)\right)
CDF	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \int_0^X \frac{1}{2xK_0(a)}\exp\left(-\frac{a}{2} \left(\frac{x}{m}+\frac{m}{x}\right)\right) \, dx
Mean	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): m\frac{K_1(a)}{K_0(a)}
Median	m
Mode	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \frac{m(\sqrt{a^2+1}-1)}{a}
Variance	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): m^2\left(1+\frac{2K_1(a)}{K_0(a)a}-\frac{K_1^2(a)}{K_0^2(a)}\right)
Skewness	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \frac{K_0^2(a)K_3(a)-3K_0(a)K_1(a)K_2(a)+2K_1^3(a)}{(K_0(a)K_2(a)-K_1^2(a))^\frac{3}{2}}
Ex. kurtosis	(see text)

In probability theory and statistics, the harmonic distribution is a continuous probability distribution. It was discovered by Étienne Halphen, who had become interested in the statistical modeling of natural events. His practical experience in data analysis motivated him to pioneer a new system of distributions that provided sufficient flexibility to fit a large variety of data sets. Halphen restricted his search to distributions whose parameters could be estimated using simple statistical approaches. Then, Halphen introduced for the first time what he called the Harmonic distribution pr Harmonic Law.. The harmonic law is a special case of the generalized inverse Gaussian distribution family when $\gamma=0$.

History

One of Halphen’s tasks,while working as statistician for Electricité de France, was the modeling of the monthly flow of water in hydroelectric stations. Halphen realized that the Pearson system of probability distributions could not be solved, it was inadequate for his purpose despite its remarkable properties. Therefore, Halphen's objective was to obtain a probability distribution with two parameters, subject a exponential decay both for large and small flows.

In 1941, Halphen decided that, in suitably scaled units, the density of X should be the same as 1/X. Taken this consideration, Halphen found the Harmonic density function. Nowadays known as an hyperbolic distribution, has been studied by Rukhin (1974) and Barndorff-Nielsen (1978).

The Harmonic Law is the only one two-parameter family of distributions that is closed under change of scale and under reciprocals, such that the maximum likelihood estimator of the population mean is the sample mean (Gauss' principle).

In 1946, Halphen realized that introducing an additional parameter, flexibility could be improved. His efforts led him to generalize the Harmonic Law to obtain the Generalized Inverse Gaussian Distribution density.

Definition

Notation

The Harmonic distribution will be denoted by Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \theta(m,a) . As a result, when a random variable X is distributed following a Harmonic Law, the parameter of scale m is the population median and a is the parameter of shape.

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): X\ \sim\operatorname{Harm}(m,a)\,

Probability density runction

The density function of the harmonic law, which depends of two parameters, has the form,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): f(x;m,a)= \frac{1}{2xK_0(a)}\exp\left(-\frac{a}{2}\left(\frac{x}{m}+\frac{m}{x}\right)\right)

where:

• Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): K_0(a)

denotes the third kind of the modified Bessel function with index  0.

• Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): m \ge 0,

• Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): a \ge 0.

Cumulative distribution function

The cumulative distribution function for the Harmonic Law does not have a closed form, and consequently it is not possible to derive an explicit expression. The cumulative distribution function must be calculated solving numerically an integral,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): F(x; m,a)=\int_0^x \frac{1}{2tK_{0}(a)}\exp\left(-\frac{a}{2}\left(\frac{t}{m}+\frac{m}{t}\right)\right) \, dt

Quantiles

The quantiles of the harmonic law are calculated with the cumulative distribution. In general does not exist a closed form except for the second quantile. We can only get the quantiles numerically.

• Q1

The first quantile, can be obtained solving the following equation expressed in terms of the integral of the probability density function:,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \int_0^{q1} \frac{1}{2xK_0(a)}\exp\left(-\frac{a}{2}\left(\frac{x}{m}+\frac{m}{x}\right)\right) \, dx= \frac{1}{4}

• Q2

In this distribution, m is the median, that is, q2 = m. Therefore,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \int_0^m \frac{1}{2xK_0(a)}\exp\left(-\frac{a}{2}\left(\frac{x}{m}+\frac{m}{x} \right)\right) \, dx= \frac{1}{2}

• Q3

Finally, the third quantile, comes from the solution of the equation,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \int_0^{q3} \frac{1}{2xK_0(a)} \exp\left(-\frac{a}{2}\left(\frac{x}{m}+\frac{m}{x}\right)\right) \, dx= \frac{3}{4}

Properties

Moments

To derive an expression for the non-central moment of order r, the integral representation of the Bessel functioncan be used.

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mu'_r = \int_0^\infty x^r f(x;m,a) \, dx= m^r \frac{K_r(a)}{K_0(a)}

where:

• r denotes the order of the moment.

Hence the mean and the succeeding three moments about it are

Order	Moment	Cumulant
1	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mu_1 = m\frac{K_1(a)}{K_0(a)}	$\mu$
2	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mu_2 = m^2\left(\frac{K_2(a)}{K_0(a)}-\frac{K_1^2(a)}{K_0^2(a)}\right)	$\sigma^2$
3	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mu_3 = m^3 \left(\frac{K_3(a)}{K_0(a)}-3\frac{K_1(a)K_2(a)}{K_0^2(a)} +2 + \frac{K_1^2(a)}{K_0^3(a)}\right)	$k_3$
4	Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \mu_4 = m^4\left(\frac{K_4(a)}{K_0(a)}-4\frac{K_1(a)K_3(a)}{K_0^2(a)} + 6\frac{K_1^2(a)K_2(a)}{K_0^3(a)} - 3\frac{K_1^4(a)}{K_0^4(a)}\right)	$k_4$

Skewness

Skewness is the third standardized moment around the mean divided by the 3/2 power of the standard deviation, we work with,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \gamma_1=\frac{\mu_3}{\mu_2^{3/2}}=\frac{K_0^2(a)K_3(a)-3K_0(a)K_1(a)K_2(a) + 2K_1^3(a)}{(K_0(a)K_2(a)-K_1^2(a))^{3/2}}

• Always $\gamma_1>0$, so the mass of the distribution is concentrated on the left.

Kurtosis

The coefficient of kurtosis is the fourth standardized moment divided by the square of the variance., for the harmonic distribution it is

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \gamma_2=\frac{\mu_4}{\mu_2^2} = \frac{K_0^3(a)K_4(a)-4K_0^2(a)K_1(a) K_3(a) + 6K_0(a) K_1^2(a) K_2(a)-3K_1^4(a)}{(K_0(a)K_2(a)-K_1^2(a))^2}

• Always Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \gamma_2>0

the distribution has a high acute peak around the mean and fatter tails.

Parameter estimation

Maximum likelihood estimation

The likelihood function is

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): L(X\mid\theta)= \prod_{i=1}^n f(x_i\mid \theta)= \prod_{i=1}^n \frac{1}{2x K_0(a)} \exp\left[-\frac{a}{2}\left(\frac{x}{m}+\frac{m}{x}\right)\right]

After that, the log-likelihood function is

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \ln L(X\mid\theta)= n\ln(2K_0(a)) + \sum_{i=1}^n \ln(x_i)+\frac{a}{2m} \sum_{i=1}^n x_i +\frac{am}{2}\sum_{i=1}^n \frac{1}{x_i}

From the log-likelihood function, the likelihood equations are,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \frac{\partial l}{\partial a} = n(\ln(2K_0(a)))'+\frac{1}{2m} \sum_{i=1}^n x_i + \frac{m}{2} \sum_{i=1}^n \frac{1}{x_i}=0

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \frac{\partial l}{\partial m} = \frac{1}{2m^2} \sum_{i=1}^n x_i +\frac{a}{2} \sum_{i=1}^n \frac{1}{x_i}= 0

To calculate the estimator $a$ only can be resolve numerically, but $m$ can see the formula.

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \hat{m}=\sqrt{\frac{\bar{H}}{\bar{H}_{-1}}}

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \sqrt{\bar{H}\bar{H}_{-1}}=\frac{K_1(\hat{a})}{K_0(\hat{a})}

Method of moments

The mean and the variance for the harmonic distribution are,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \begin{cases} \mu = m\frac{K_1(a)}{K_0(a)} \\ \sigma^2 = m^2(1+\frac{2K_1(a)}{K_0(a)a}-\frac{K_1^2(a)}{K_0^2(a)}) \end{cases}

Note that,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \sigma^2 = \mu^2(\frac{2K_0(a)}{K_1(a)})^2+\frac{2K_0(a)\mu^2}{K_1(a)a}- \mu^2

The method of moments consists in to solve the following equations:

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \begin{cases} \bar{H}=m\frac{K_1(a)}{K_0(a)} \\ s^2= \bar{H}^2(\frac{2K_0(a)}{K_1(a)})^2+\frac{2K_0(a)\bar{H}^2}{K_1(a)a}- \bar{H}^2 \end{cases}

where $s^2$ is the sample variance and $\bar{H}$ is the sample mean. Solving the second equation we obtain $\hat{a}$, and then we calculate $\hat{m}$ using,

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): \hat{m}=\frac{\bar{H}K_0(\hat{a})}{K_1(\hat{a})}

Related distributions

The harmonic law is a sub-family of the generalized inverse Gaussian distribution. The density of GIG family have the form

Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): f(x\mid m,\gamma)= \frac{x^{\gamma-1}}{2m^\gamma K_\gamma(a)}\exp\left[-\frac{a}{2} \left(\frac{x}{m}+\frac{m}{x}\right)\right]

The density of the Generalized Inverse Gaussian Distribution family corresponds to the Harmonic Law when $\gamma=0$.

When $a$ tends to infinity, the harmonic law can be approximated by a normal distribution. This is indicated by demonstrating that if $a$ tends to infinity, then Failed to parse (Missing <code>texvc</code> executable. Please see math/README to configure.): U=\sqrt{a}\left(\frac{X}{m}-1\right) , which is a linear transformation of X, tends to a normal distribution ($N(0,1)$).

This explains why the normal distribution can be used successfully for certain data sets of ratios.

Another related distribution is the log-harmonic law, which is the probability distribution of a random variable whose logarithm follows an harmonic law.

This family has an interesting property, the Pitman estimator of the location parameter does not depend on the choice of the loss function. Only two statistical models satisfy this property: One is the normal family of distributions and the other one is a three-parameter statistical model which contains the log-harmonic law.

See also

• Normal distribution
• Generalized inverse Gaussian distribution