Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. It's 100% free, no registration required.

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

How do I fit the parameters of a t-distribution, i.e. the parameters corresponding to the 'mean' and 'standard deviation' of a normal distribution. I assume they are called 'mean' and 'scaling/degrees of freedom' for a t-distribution?

The following code often results in 'optimization failed' errors.

library(MASS)
fitdistr(x, "t")

Do I have to scale x first or convert into probabilities? How best to do that?

Thanks!

share|improve this question
1  
It fails not because you have to scale parameters, but because optimizer fails. See my answer below. – Sergey Bushmanov Feb 10 at 15:26
up vote 4 down vote accepted

fitdistr uses maximum-likelihood and optimization techniques to find parameters of a given distribution. Sometimes, especially for t-distribution, as @user12719 noticed, the optimization in the form:

fitdistr(x, "t")

fails with an error.

In this case you should give optimizer a hand by providing starting point and lower bound to start searching for optimal parameters:

fitdistr(x, "t", start = list(m=mean(x),s=sd(x), df=3), lower=c(-1, 0.001,1))

Note, df=3 is your best guess at what an "optimal" df could be. After providing this additional info your error will be gone.

Couple of excerpts to help you better understand the inner mechanics of fitdistr:

For the Normal, log-Normal, geometric, exponential and Poisson distributions the closed-form MLEs (and exact standard errors) are used, and start should not be supplied.

...

For the following named distributions, reasonable starting values will be computed if start is omitted or only partially specified: "cauchy", "gamma", "logistic", "negative binomial" (parametrized by mu and size), "t" and "weibull". Note that these starting values may not be good enough if the fit is poor: in particular they are not resistant to outliers unless the fitted distribution is long-tailed.

share|improve this answer
    
Both answers (Flom and Bushmanov) are helpful. I'm picking this one, because it makes it more explicit that with the right initial values and constraints 'fitdistr' optimization converges. – user12719 Feb 10 at 23:01

In the help for fitdistr is this example:

fitdistr(x2, "t", df = 9)

indicating that you just need a value for df. But that assumes standardization.

For more control, they also show

mydt <- function(x, m, s, df) dt((x-m)/s, df)/s
fitdistr(x2, mydt, list(m = 0, s = 1), df = 9, lower = c(-Inf, 0))

where the parameters would be m = mean, s = standard deviation, df = degrees of freedom

share|improve this answer
    
I guess I am confused about the parameters of a t-distribution. Does it have 2 (mean,df) or 3 (mean, standard deviation, df) parameters? I was wondering if one could fit the parameter 'df'. – user12719 Dec 12 '12 at 21:10
1  
@user12719 The Student's-t distribution has three parameters: location, scale and degrees of freedom. They are not referred as mean, standard deviation and df because the mean and the variance of this distribution depend on the three parameters. Also, they do not exists in some cases. Peter Flom is fixing the df but this can be considered as an unknown parameter as well. – user10525 Dec 12 '12 at 21:25
1  
@PeterFlom In the case of the Cauchy distribution it is explicit that m and s are the location and scale. I agree the notation m and s suggests that they represent the mean and standard deviation, respectively. But this may just be a simplification of \mu and \sigma as well. +1 long ago, by the way. – user10525 Dec 12 '12 at 22:09
1  
@PeterFlom Does this citation from the R's help file imply that df is always 9 for students distribution? Don't you think df should be estimated as well? Actually, the absence of df is the cause for the error, and the right answer should provide some recipe to finding it. – Sergey Bushmanov Feb 9 at 17:23
1  
@PeterFlom BTW, if you read the help file couple of lines above your citation, you'll find why df=9 is good in their example and irrelevant here. – Sergey Bushmanov Feb 9 at 17:30

MASS, the book (4th edition, page 110) advices against trying to estimate $\nu$, the degrees of freedom parameter in the t distribution with maximum likelihood (with some literature references I will add in here when I have the book with me).

The reason is that the likelihood function for $\nu$ based on the t density function, may be unbounded! and will in those cases not give a well defined maximum. Let us see at an artificial example where location and scale is known (as the standard $t$-distribution) and only the degrees of freedom is unknown. Below is some R code, simulating some data, defining the loglikelihood function and plotting it:

set.seed(1234)
n <- 10
x <- rt(n,  df=2.5)

make_loglik  <-  function(x)
    Vectorize( function(nu) sum(dt(x, df=nu,  log=TRUE)) )

loglik  <-  make_loglik(x)
plot(loglik,  from=1,  to=100,  main="loglikelihood function for df     parameter", xlab="degrees of freedom")
abline(v=2.5,  col="red2")

enter image description here

If you play around with this code, you can find some cases where there is a well-defined maximum, especially with the sample size $n$ large. But is the maximum likelihood estimator then any good?

Let us try some simulations:

t_nu_mle  <-  function(x) {
    loglik  <-  make_loglik(x)
    res  <-  optimize(loglik, interval=c(0.01, 200), maximum=TRUE)$maximum
    res   
}

nus  <-  replicate(1000, {x <- rt(10, df=2.5)
    t_nu_mle(x) }, simplify=TRUE)

> mean(nus)
[1] 45.20767
> sd(nus)
[1] 78.77813

Showing the estimation is very unstable (looking at the histogram, a sizeable portion of the estimated values is at the upper limit given to optimize of 200).

Repeating with a larger sample size:

nus  <-  replicate(1000, {x <- rt(50, df=2.5)
    t_nu_mle(x) }, simplify=TRUE)
> mean(nus)
[1] 4.342724
> sd(nus)
[1] 14.40137

which is much better, but the mean is stil way above the true value of 2.5.

Then remember that this is a simplified version of the real problem when location and scale parameters also have to be estimated.

If the reason of using the $t$-distribution is to "robustify", then estimating $\nu$ from the data well may destroy the robustness.

share|improve this answer
1  
Your conclusion that the problems of estimating df might actually work against the reason to choose a t-distribution in the first place (i.e. robustness) is thought provoking. – user12719 Feb 12 at 20:00

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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