0
$\begingroup$

I have two questions, a conceptual one and a practical one (that are closely related). And just as a note, I'm not super familiar with this level of stats (much more comfortable with simple linear regression and ANOVA), so apologies if there are oversimplifications/obvious conceptual issues at play here.

First: There are clearly issues with using means to analyze reaction time data because RT data are not normally distributed. An appropriate distribution that's frequently referenced is the ex-Gaussian distribution. Because I wanted to implement an analysis using an appropriate distribution without transforming the data I looked into using the GAMLSS package to run the analysis. But then I go back to the literature and what I see is that the ex-Gaussian parameters (mu, sigma, and tau) are being extracted from the gam model and then being used in a GLMM. So the question is: why not just report and interpret the gam model with an ex-Gaussian distribution? Is there anything incorrect about doing so in order to examine these data with a non-normal distribution?

Second: a brief rundown of the study then the questions (sample output is at the end). Participants completed a simple RT task at two timepoints at two different locations. The RT task has a condition called valence with three levels. Participants are randomized to watch 1 of 2 videos at both timepoints. Participants complete a measure called STAI once. Variables: location (w/in factor with two levels); valence (w/in factor with three levels); video (between factor with two levels); STAI (continuous variable). Random intercepts are included for individual participants.

There are a few things I don't understand about GAMLSS: 1) The output is different when I include different estimation parameters (so if I include only mu in the formula e.g. gamlss(formula = RT~location, data=data) it's different from gamlss(RT~location, sigma.fo = ~location, nu.fo = ~location, data=data) Any thoughts on why this is the case? Are fit statistics reliable between these two models? I post an example output with the mu and nu formulas only at the end. 2) This is a continuation of the first part: my understanding of interpreting the parameters is that they each refer to a distinct component (mu refers to mean; sigma to sd; tau to the exponential component of the mean and sd). So if location (which strongly influences RT in every model I try and run except when I included the tau component in the analysis) is significant for tau, but not mu - does this mean that the real explanation for the difference of location on RTs is not global (i.e. average) but somehow explained by the exponential convolution of the mean and sd? How would I explain this?

Below I provide an example output of the model I ran first with all formulas, then only the mu and nu formula (gamlss refers to nu, but it's actually equivalent to tau in the context of the exGAUS distribution, I think?).

Family:  c("exGAUS", "ex-Gaussian") 

Call:  
gamlss(formula = RT ~ STAI_T_tot + location + video + valence +  
    random(subID), sigma.formula = ~STAI_T_tot + location + video +  
    valence + random(subID), nu.formula = ~STAI_T_tot + location +  
    video + valence + random(subID), family = exGAUS(), data = na.omit(data)) 


Fitting method: RS() 

-------------------------------------------------------------------
Mu link function:  identity
Mu Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)    606.19898    2.48220 244.219  < 2e-16 ***
STAI_T_tot      -0.78644    0.05277 -14.905  < 2e-16 ***
location1        1.03000    1.11498   0.924   0.3556    
video1           8.08323    1.11894   7.224 5.09e-13 ***
valenceneutral   3.37995    1.35687   2.491   0.0127 *  
valencepos       2.82100    1.35826   2.077   0.0378 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
Sigma link function:  log
Sigma Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)     4.6897578  0.0293722 159.667  < 2e-16 ***
STAI_T_tot     -0.0015505  0.0006102  -2.541  0.01106 *  
location1       0.0236655  0.0131948   1.794  0.07289 .  
video1          0.0384450  0.0131913   2.914  0.00356 ** 
valenceneutral  0.0185930  0.0161728   1.150  0.25029    
valencepos      0.0255943  0.0161156   1.588  0.11225    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
Nu link function:  log 
Nu Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)     5.6664104  0.0150130 377.432   <2e-16 ***
STAI_T_tot     -0.0030183  0.0003188  -9.467   <2e-16 ***
location1       0.1074412  0.0070122  15.322   <2e-16 ***
video1          0.0640823  0.0070783   9.053   <2e-16 ***
valenceneutral -0.0107050  0.0085951  -1.245    0.213    
valencepos     -0.0102235  0.0085242  -1.199    0.230    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
NOTE: Additive smoothing terms exist in the formulas: 
 i) Std. Error for smoothers are for the linear effect only. 
ii) Std. Error for the linear terms may not be reliable. 
-------------------------------------------------------------------
No. of observations in the fit:  80795 
Degrees of Freedom for the fit:  436.9502
      Residual Deg. of Freedom:  80358.05 
                      at cycle:  20 

Global Deviance:     1115911 
            AIC:     1116784 
            SBC:     1120848 

`Family:  c("exGAUS", "ex-Gaussian") 

Call:  
gamlss(formula = RT ~ STAI_T_tot + location + video + valence +  
    random(subID), nu.formula = ~STAI_T_tot + location + video +  
    valence + random(subID), family = exGAUS(), data = na.omit(data)) 


Fitting method: RS() 

-------------------------------------------------------------------
Mu link function:  identity
Mu Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)    540.62127    2.58243 209.346  < 2e-16 ***
STAI_T_tot      -0.15510    0.05369  -2.889  0.00387 ** 
location1       -1.34292    1.17034  -1.147  0.25119    
video1         -15.72136    1.17199 -13.414  < 2e-16 ***
valenceneutral   1.37373    1.43494   0.957  0.33840    
valencepos       0.60667    1.43042   0.424  0.67148    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
Sigma link function:  log
Sigma Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 4.418776   0.006184   714.6   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
Nu link function:  log 
Nu Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)     5.8552886  0.0152149 384.838   <2e-16 ***
STAI_T_tot     -0.0045268  0.0003249 -13.934   <2e-16 ***
location1       0.1100487  0.0071404  15.412   <2e-16 ***
video1          0.0889630  0.0072172  12.327   <2e-16 ***
valenceneutral -0.0034509  0.0087375  -0.395    0.693    
valencepos     -0.0029106  0.0086685  -0.336    0.737    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

-------------------------------------------------------------------
NOTE: Additive smoothing terms exist in the formulas: 
 i) Std. Error for smoothers are for the linear effect only. 
ii) Std. Error for the linear terms may not be reliable. 
-------------------------------------------------------------------
No. of observations in the fit:  80795 
Degrees of Freedom for the fit:  294.7012
      Residual Deg. of Freedom:  80500.3 
                      at cycle:  20 

Global Deviance:     1119469 
            AIC:     1120059 
            SBC:     1122799 `
$\endgroup$

1 Answer 1

2
$\begingroup$

There is a paper covering the package, but I'm not sure how much help to you it might be.

Based on my understanding the purpose of GAMLSS is to be able to run different models of distribution parameters in the same global model. This is similar to mixed-effects regression where the model is hierarchical, i.e., that equations for parameters at different levels are fit in one model in order to better model sources of error (e.g., dependency between repeated measures responses). But but GAMLSS you are modeling specific explanatory effects on specific parameters from whatever distribution you have chosen (in your case, the ex-Gaussian).

That works well because in the RT literature using the ex-Gaussian as the model of cognition, there are theoretical reasons to believe that experimental manipulations will affect one or more parameters by not the others (but see Matzke & Wagenmakers, 2009). In your case, one might theorize that, say, valence only affects $\mu$, but not the other parameters (I am just picking at random here). So then you'd specify an equation only for $\mu$ but not the other variables. Everything else would go in the first formula call, as you have done.

If you do not enter an explicit equation for a parameter, the function just fits an "intercept", i.e., an estimate of that parameter without any knowledge of factors or covariates.

[1]: Matzke, D., & Wagenmakers, E.-J. (2009). Psychological interpretation of the ex-Gaussian and shifted Wald parameters: A diffusion model analysis. Psychonomic Bulletin & Review, 16(5), 798–817. http://doi.org/10.3758/PBR.16.5.798

$\endgroup$
1
  • $\begingroup$ Also, to put sigma and nu back on the scale of seconds, you must exponentiate the estimate for the intercept, as the link function is log. $\endgroup$
    – Rick Hass
    Commented Jun 23, 2017 at 20:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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