3
$\begingroup$

Let $x$ be a $\alpha-$stable distributed random variable of parameters $\alpha,\beta,c,\mu$. When $\alpha \gt 1$ I can estimate the location parameter $\mu$ of the distribution as

$\mu=E[x]$

But how do estimate $\mu$ when $0\lt \alpha \le 1$ ?

$\endgroup$
3
  • 1
    $\begingroup$ "$E[x]$", if it even exists, is an unknown theoretical quantity, not an estimator. I believe you intended that to be the sample mean. $\endgroup$
    – whuber
    May 14 '15 at 18:51
  • $\begingroup$ yes, with E[x] I intend the sample mean. $\endgroup$
    – emanuele
    May 14 '15 at 19:51
  • $\begingroup$ Here is a paper on arXiv: arxiv.org/pdf/1706.09756.pdf $\endgroup$ Nov 10 '18 at 17:32
4
$\begingroup$

Estimating $\mu$ using the sample mean is only a reasonable option if $\beta \equiv 0$, i.e., a symmetrical distribution. Otherwise the sample average, $\bar{x}$, is severely biased for $\mu$, the location parameter -- not to be confused with the expected value of a stable distribution, $E(X)$, which is $\neq \mu$ if $\beta \neq 0$ and $\alpha < 2$ (for $\alpha = 2$, distribution is always a Normal distribution independent of $\beta$).

For what follows I assume $\beta = 0$. For skewed distributions your best bet is probably a likelihood based optimization / posterior calculation.

For $\alpha < 2$ a natural candidate to estimate the location in a more robust manner than the sample mean, is to use the sample median. Another alternative is to use heavy-tail Lambert W x Normal distribution estimates and use $\hat{\mu}$ from the Lambert W x Normal fit as an estimate of the location of your data (from the stable distribution). See an example of this here, which compares estimators for the location of a Cauchy ($\alpha = 1$) -- including sample mean, sample median, and the Lambert W x Normal options.

To address your specific question about $\alpha < 1$ I replicated the simulation from that post, and simulated samples from a stable distribution with $\alpha = 0.75$ and $\mu = 0$. Results show that in this case Lambert W x Normal estimators even outperform the sample median.

library(alphastable)
library(LambertW)

LocationEstimators <- function(x.sample) {
  out <- 
    c(mean = mean(x.sample),
      median = median(x.sample))

  # Lambert W x Gaussian estimates for heavy tails ('h')
  igmm.tau <- LambertW::IGMM(x.sample, "h")$tau
  beta.hat <- igmm.tau[1:2]
  names(beta.hat) <- c("mu", "sigma")
  mle.lambertw <- LambertW::MLE_LambertW(x.sample, distname = "normal", 
                               theta.init = LambertW::tau2theta(igmm.tau, 
                                                      beta = beta.hat),
                               type = "h",
                               return.estimate.only = TRUE)
  out <- c(out, igmm.tau["mu_x"], mle.lambertw["mu"])
  names(out)[3:4] <- c("igmm.LambertW", "mle.LambertW")
  return(out)
}

sim_est <- function(n, ...) {
  yy <- alphastable::urstab(n=n, param=0, ...)
  return(LocationEstimators(yy))
}

These functions can now be used to run the simulation study with $n = 100$ replications and a sample size of 1000 each.

# Simulate and look at bias, std dev, and MSE
nsim <- 100
num.samples <- 1000
set.seed(nsim)
est = t(replicate(n=nsim, sim_est(n=num.samples, alpha=0.75, beta=0, sigma=1., mu=0.)))
colMeans(est)
         mean        median igmm.LambertW  mle.LambertW 
-1.147466e+02 -5.030865e-03 -2.193856e-03 -2.580000e-03 

apply(est, 2, sd)
         mean        median igmm.LambertW  mle.LambertW 
 764.12848769    0.04120362    0.07423480    0.03601814 

# RMSE (since true location = 0)
sqrt(colMeans(est^2))
         mean        median igmm.LambertW  mle.LambertW 
 768.90845365    0.04130460    0.07389526    0.03593034 

As expected, sample mean is useless; median and Lambert W estimators all provide (effectively) unbiased estimates of $\mu = 0$, with Lambert W MLE < median < Lambert W IGMM in terms of standard errors / MSE.

A violin plot of these 3 estimators makes this clear as well:

library(ggplot2)
library(reshape2)
theme_set(theme_bw(18))
est.m <- melt(est)
colnames(est.m) <- c("sim.id", "estimator", "value")

# remove 'mean' for good scaling in plots
est.m <- subset(est.m, estimator != 'mean')
ggplot(est.m,
       aes(estimator, value, fill = estimator)) +
  geom_violin() +
  geom_hline(yintercept = 0, size = 1, linetype = "dashed", 
             colour = "blue") +
  theme(legend.position = "none",
        axis.text.x = element_text(angle = 90)) + ylab("location estimate")

location_estimator_comparison

$\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.