1
$\begingroup$

Let's say we have a log-logistic random variable $X$ with probability density function:

$$f(x)=\frac{(\beta/\alpha)(x/\alpha)^{\beta-1}}{(1+(x/\alpha)^{\beta})^{2}}$$

where $\alpha>0$, $\beta>0$ and support on $x\in[0,\infty)$.

This random variable has mean and variance defined as follows:

$$\mathbb{E}[X]=\frac{\alpha(\pi/\beta)}{\text{sin}(\pi/\beta)}$$ $$\text{Var}[X]=\alpha^{2}\bigg[\frac{2(\pi/\beta)}{\text{sin}(2(\pi/\beta))}-\frac{(\pi/\beta)^{2}}{\text{sin}^{2}(\pi/\beta)}\bigg]$$

or, more generally, with $k$th raw moment defined for $k<\beta$: $$\mathbb{E}[X^{k}]=\alpha^{k}\frac{k(\pi/\beta)}{\text{sin}(k\pi/\beta)}$$

Now, is there an clever way to estimate the parameters $(\alpha,\beta)$ such that it provides a desired mean and variance?

As far as I can tell, we can arrive at an equation in terms of either parameter, but solving requires some form of optimization. I was curious to know if there are any methods to avoid such optimization?

$\endgroup$
  • $\begingroup$ Are you familiar with method of moments estimators? $\endgroup$ – AdamO Oct 30 '17 at 4:54
  • $\begingroup$ @AdamO yes, I couldn't find any MOM estimators for the log-logistic $\endgroup$ – Ed P Oct 30 '17 at 4:55
  • 2
    $\begingroup$ Variance only exists for $\beta>2$, and in that case $E(X)^2/E(X^2)$ s a monotonic increasing function only of $\beta$. $\endgroup$ – Glen_b -Reinstate Monica Oct 30 '17 at 7:23
  • 3
    $\begingroup$ Hmm, maybe that's the hard way around. Let $\theta=\pi/\beta$. Then $(\text{Var}(X)/E(X)^2) + 1=E(X^2)/E(X)^2 = \tan(\theta)/\theta$ for $0<\theta<\pi/2$ is monotonic in $\theta$. This is relatively easy to solve numerically and then take back to an estimate of $\beta$, and then solve for $\alpha$. $\endgroup$ – Glen_b -Reinstate Monica Oct 30 '17 at 8:09
1
$\begingroup$

Note that $E(X^k)$ has $\alpha$ present only as a multiplicative term in $\alpha^k$.

As a result $E(X^a)^b/[E(X^b)^a]$ is a function of $\beta$ alone. Since $E(X^k)$ only exists for $\beta>k$ it's best to use lowest moments.

In particular, $E(X^2)/(E(X)^2)$ or $\operatorname{Var}(X)/(E(X)^2)$ or convenient functions of them may be used to solve for $\beta$ (at least numerically). Once $\beta$ is known, it is simple to obtain $\alpha$ from $E(X)$.

Now let's look at some convenient ways to obtain $\beta$.

Note that $\frac{2\sin^2(x)}{\sin(2 x)} = \frac{2\sin^2(x)}{2\sin(x)\cos(x)}= \tan(x)$ (on $0<x<\pi/2$).

Let $\theta=\pi/\beta$. Then $E(X)^2/E(X^2)=\theta/\tan(\theta)$ for $0<\theta<\pi/2$ is monotonic decreasing in $\theta$, and has the added benefit that the function is bounded in this range. This should be reasonably simple to solve for $\theta$ using standard numerical root finding methods. Because $E(X)^2/E(X^2)$ will always lie between $0$ and $1$, this should always yield a possible solution (one with $\beta>2$)

plot of function x/tan(x) on 0 to \pi/2; flat (at 1) at x=0 then more rapidly decreasing as it gets toward x=\pi/2 - seems close to quadratic

However, we can get very close to to the solution -- and have bounds on it -- with just a little more manipulation.

Note that $\sqrt{1-\frac{x}{\tan(x)}}$ is very close to linear:

plot of function \sqrt{1-\frac{x}{\tan(x)}} or 0<x<pi/2

In that interval it is bounded between the lines $x/\sqrt{3}$ and $\frac{2x}{\pi}$ -- which slopes only differ by about 10%.

So the steps would be:

  • Compute $k=\sqrt{\frac{\sigma^2}{\sigma^2+\mu^2}}$.

  • Solve the equation $\sqrt{1-\frac{\theta}{\tan(\theta)}}-k=0$ for $\theta$ via a convenient root-finding algorithm, with initial bounds $\frac{\pi}{2}k<\theta<\sqrt{3}k$. This should be very rapid.

  • Obtain $\beta=\pi/\theta$

  • Calculate $\alpha=\mu\frac{\sin\theta}{\theta}$.

Example: Let's say we want to figure out the parameters of a log-logistic distribution that would have $\mu=10$ and $\sigma=10$.

Then $k=\sqrt{\frac{10^2}{10^2+10^2}}=\sqrt{2}$. Using R:

f=function(x,k) {sqrt(1-x/tan(x))-k}
mu=10
sigma=10
k=sqrt(sigma^2/(sigma^2+mu^2))
f0 = uniroot(f,lower=k*pi/2,upper=min(pi/2,k*sqrt(3)),k=k)
f0
$root
[1] 1.165561

$f.root
[1] -3.517283e-08

$iter
[1] 3

$init.it
[1] NA

$estim.prec
[1] 6.103516e-05

So $\theta\approx 1.165561$.

 theta=f0$root
 beta=pi/theta
 alpha=10*sin(theta)/theta
 alpha;beta
[1] 7.884698
[1] 2.695348

Now let's check we actually got the mean and sd we were aiming for:

 alpha*pi/beta/sin(pi/beta)
[1] 10
 alpha*sqrt(2*pi/beta/sin(2*pi/beta)-(pi/beta)^2/sin(pi/beta)^2)
[1] 9.999999

About 7 digits of accuracy that time (we were probably a bit lucky getting the root quite so close to 0 that time; the estimated precision suggests it should not have worked quite that well).

Looks good; it seems to be doing what it should.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.