I have several response times from participants in a mental rotation task. I'm really puzzled how to get the best distribution model to fit my data. I had checked several ones with the fitdistrplus::fitdist function and have two good candidates (log-normal & gamma).
I would like to check for an ex-gaussian distribution fit (GAMLSS package), but the fitdist function ask for "start" parameters and I have no idea how to get them: any help would be much appreciated.
I did a quick google and found the mesgauss function in the retimes package.
I give a simple example of estimatinat ex-gaussian parameters from sample data.
library(retimes)
rtdata <- c(.5, .7, .3, 1, 2, .5, .55)
retimes::mexgauss(rtdata)
This returns:
mu sigma tau
0.2838158 0.2668972 0.5090413
Regarding starting values, typically these just need to be ball-park estimates of what the parameters might be.
To answer the question relating to start values for the parameters for use with fitdist:
I would like to check for an ex-gaussian distribution fit (GAMLSS package), but the fitdist function ask for "start" parameters and I have no idea how to get them: any help would be much appreciated.
Where possible, knowledge of the behavior of similar reaction times can be used. This might be particularly helpful when the sample skewness (see below) is not suitable for obtaining starting values, since skewness in particular might be anticipated to reasonably similar across different but related circumstances.
As Jeromy says, start values generally just need to be rough ballpark estimates, so if you have some experience of fits you may be able to supply suitable guesses. But if you don't have any suitable idea, there's still something you can do.
The Wikipedia page on the ExGaussian distribution offers some possibilities.
Firstly, it gives method-of-moments-related estimators based on the sample mean, variance and skewness (labelled $m, s^2$, and $\gamma_1$ there):
$${\displaystyle {\hat {\mu }}=m-s\left({\frac {\gamma _{1}}{2}}\right)^{1/3},}$$
$${\displaystyle {\hat {\sigma }^{2}}=s^{2}\left[1-\left({\frac {\gamma _{1}}{2}}\right)^{2/3}\right],}$$
$${\displaystyle {\hat {\tau }}=s\left({\frac {\gamma _{1}}{2}}\right)^{1/3}.}$$
The main thing to note there is that the $\tau$ (scale) parameter there is the same as the $\nu$ parameter in the GAMLSS function exGAUS (and is the inverse of the $\lambda$ rate-parameter in the definition of the distribution at the top of the Wikipedia article). The easy way to do this is by starting with $\hat{\tau}$:
$${\displaystyle {\hat {\tau }}=s\left({\frac {\gamma _{1}}{2}}\right)^{1/3}.}$$
then
$$ \hat {\mu }=m-\hat{\tau}$$
$$\hat {\sigma}^{2}=s^{2}-\hat{\tau}^2\,.$$
However, you need to be careful that you don't end up with implausible values for the parameters which might happen if the sample skewness were negative (leading to an impossible value for $\tau$) or if the sample skewness were very large (perhaps leading to issues with the other two parameters, such as a negative variance) -- this must be checked.
On the other hand problems here should also be seen as a warning; if those things happen with the method of moments, there's strong potential for the likelihood function to not behave nicely as well.
The article also mentions the range of the nonparametric skew (a skewness measure based on the second Pearson skewness), which relates the mean, median and standard deviation. This might be used in place of the moment-skewness (with some suitable modification) to obtain an estimate of the scale parameter of the exponential component.
An alternative would be to look at quantile-matching (such as a calculation based on a high, low and middling quantile) or at estimates from quantile-based measures of location, scale and skewness, but the quantiles are not in general simple functions of the parameters; the existence of quantile and cdf functions in GAMLSS make this potentially viable as a strategy, however. This would be more complex but more likely to lead to feasible values of the parameters, at least in the case where the skewness is very large.
The notion of finding a "best" distribution -- and in particular searching for one with a given set of data -- is not especially useful. No simple distribution will be a perfect fit (no simple model will be exactly correct). Models are useful representations rather than perfect descriptions. If one has a simple model which usually fits fairly well and has some theoretical basis by which its parameters and relationships between distributions based on different estimates can be understood, that's a good reason to use it. We should not be overly concerned about a "best" model for a particular data set so much as an adequate model across a variety of data of similar kinds.
In relation to the title -- fitdist (from package fitdistrplus, or fitdistr in MASS) will give parameter estimates for a fitted distribution but doesn't really tell you whether the fit is good enough for some particular purpose or other (i.e. it doesn't really check the quality of the fit in a practical sense).
