# Transformation to fit gamma distribution for glm

The data simulated below has a maximum value of 4 and is interestingly skewed. The maximum of 4 is a limitation imposed by the instrument used and the data is semi-discrete, i.e., there are a reasonably large number of numbers it could be between -4 and 4. Because of the shape of the data, I thought about transforming it so it would approximate a gamma distribution:

It is limited to this range in this instance because it is a signal detection measure (d prime http://en.wikipedia.org/wiki/D%27) and the accuracy we have for this particular measure limits us to +-4. It is skewed like this because one population does not very often get false positives and will generally get more hits while the other populations often do get false positives and less hits.

set.seed(69)
g1<-rnorm(700,0,1); g2<-rnorm(100,-0.5,1.5); g3<-rnorm(100,-1,2.5)
gt<-data.frame(score=c(g1, g2, g3), fac1=factor(rep(c("a", "b", "c"), c(700, 100, 100))), fac2=ordered(rep(c(0,1,2), c(3,13,4))))
gt$score<-with(gt, ifelse(fac2 == 0, score, score-rnorm(1, 0.5, 2))) gt$score<-with(gt, ifelse(fac2 == 2, score-rnorm(1, 0.5, 2), score))
gt$score<-round(with(gt, ifelse(score>0, score*-1, score)), 1)+4 gt$score<-with(gt, ifelse(score < -4, -4, score))
gt$cov1<-with(gt, score + rnorm(900, sd=40))/40 hist(gt$score)
gt$score2<-with(gt, 4-score+0.0000001) #Gamma distribution can't have 0s (and is positive skewed???) hist(gt$score2)

glm1<-glm(score2~cov1+fac1*fac2, family="Gamma", data=gt)


This is quite new territory for me.
1. Is this a reasonable thing to do?
2. Are there other distributions I might try and compare (exponential perhaps)?

Update:
After some comments below, I investigated beta regression using the betareg package in R. It gave me skewed residuals:

gt$scorer<-with(gt, (score--4)/(4--4)) gt$scorer<-with(gt, (scorer*(length(scorer)-1)+0.5)/length(scorer))
b1 <- betareg(scorer ~ cov1 + fac1 * fac2, data=gt)
plot(density(resid(b1))) #Strange residuals, even straight lm looks better


So I had a look at a quasibinomial regression and it gave me smaller and better looking residuals:

glm2 <- glm(scorer~cov1 + fac1 * fac2, data=gt, family="quasibinomial")
plot(density(resid(g1))) #Better residuals


Are the residuals good enough to go on in this case?
Or is the fact that d', while based upon T/F, is not a binary variable, a serious issue?

Edit 3: d' clarification The below is an example of my d' scores, with the rough distributional qualities and similar raw scores for hits and false positives.

hitrate<-sample(0:16, 100, replace=T, prob=c(rep(0.02,11), 0.025, 0.05, 0.1, 0.2, 0.3, 0.2))/16
hitrate<-ifelse(hitrate==1, 31/32,hitrate); hitrate<-ifelse(hitrate==0, 1/32,hitrate)
farate<-sample(0:32,100, replace=T, prob=c(0.7,0.1,0.05,0.05,0.05,0.02,rep(0.001, 27)))/32
farate<-ifelse(farate==0, 1/64,farate); farate<-ifelse(farate==1, 63/64,farate)

dprime<-round(qnorm(hitrate) - qnorm(farate),1)
plot(density(dprime))

• Why would you think your data should follow the gamma distribution, which has an infinite support (all $x>0$)? To find the right distribution, tell us about the data generating process: what makes it restricted to the $(-4,4)$ range, and why is it skewed? – Aniko Nov 7 '11 at 20:14
• Updated question to include more info on the data. It appears the gamma distribution is certainly not correct for this instance. – Matt Albrecht Nov 8 '11 at 5:06
• This is a very interesting problem, but I need further clarification. The following is based on the signal detection link you included. Is this correct? You send $n_T$ true signals, and get $k_T$ detections for a hit rate of $h=k_T/n_T$; then you send false signals for a false alarm rate of $f=k_F/n_F$, and finally you get $d'=\Phi^{-1}(h) - \Phi^{-1}(f)$, which is the score you are analyzing. If this is correct, would you have any of the intermediate values available? – Aniko Nov 9 '11 at 15:01
• Yes, its the same signal detection calculation. I've updated the question to show the rough process; from raw scores to final d' scores. Is this what you were after? – Matt Albrecht Nov 9 '11 at 15:53

A gamma distribution definitely doesn't make sense for you data. Gamma takes support on the entire real line and is always skewed right. The example data you provide in your code would be horrible data to try to fit a gamma to.

It would definitely be nicer to know more about the data generation process. But one thing that comes to mind is you could scale and shift the data to be constrained between 0 and 1 and then attempt to use a Beta distribution to model that. Once again it would be better to know more about your data but a Beta is one of the few well known parametric distributions that is bounded below and above.

However, it seems you want to do some sort of a regression. Have you tried to fit the regression assuming a normal error term and examining the residuals? A lot of people assume that the data itself needs to be normally distributed for a linear regression to work but typically we place the assumption on the error term and depending on the values your covariates take this can lead to a skewed distribution for your response variable.

• Thanks. So gamma distributions are really only for situations where the possibility could be any number, I made a similar error previously with respect to negative binomials. A beta binomial type distribution as you suggest sounds good. The residuals show a slightly skewed distribution in this case. However, sometimes the measures I'm interested in are not so well behaved for the normal regression. I have updated the question to give more info about the data. – Matt Albrecht Nov 8 '11 at 5:04
• Oops, realised you meant the beta distribution not the beta-binomial. – Matt Albrecht Nov 8 '11 at 5:25
• I wouldn't rely heavily on using considerations of support to select a probability model. After all, normal distributions--which, like Gammas, have infinite support--are frequently used, with great success, to model distributions that must have finite support (such as the sampling distribution of binary responses in a survey question). – whuber Nov 8 '11 at 14:49
• Very true. Which is why I wondered what happened with the residuals for a normal linear regression. Using their example data it didn't seem like a gamma was very appropriate. – Dason Nov 8 '11 at 17:28

I was going to suggest using binomial regression with a probit link on the original hit/false alarm counts, but before writing out the details I googled for this idea. Apparently, someone else thought of it already (there goes a publication :)

Here is the reference (in case the link is restricted), there is also an R package sensR to go with it:

PB Brockhoff, RH Bojesen Christensen (2010) Thurstonian models for sensory discrimination tests as generalized linear models, Food Quality and Preference, 21(3), 330-338.

• Thanks, will probably be back with a question later after reading the ref. Would give an extra +1 for the reference if I could. – Matt Albrecht Nov 10 '11 at 1:44