27

My dependent variable shown below doesn't fit any stock distribution that I know of. Linear regression produces somewhat non-normal, right-skewed residuals that relate to predicted Y in an odd way (2nd plot). Any suggestions for transformations or other ways to obtain most valid results and best predictive accuracy? If possible I'd like to avoid clumsy categorizing into, say, 5 values (e.g., 0, lo%, med%, hi%, 1).

enter image description here

enter image description here

amoeba
  • 104,745
rolando2
  • 12,511
  • 7
    You would be better off telling us about these data and where they came from: something has clamped a distribution that naturally extends beyond the $[0,1]$ interval. It's possible you have used some measurement method or statistical procedure that's not quite appropriate for your data. Trying to patch up such a mistake with sophisticated distribution-fitting techniques, nonlinear re-expressions, binning, etc., would just compound the error, so it would be nice to circumvent the problem altogether. – whuber Feb 06 '13 at 21:10
  • 3
    @whuber - A good thought, but the variable was created through a complex bureacratic system which is unfortunately set in stone. I'm not at liberty to disclose the nature of the variables involved here. – rolando2 Feb 06 '13 at 21:30
  • 2
    Okay, it was worth a shot. I'm thinking that instead of transforming the data, you might still want to recognize the clamping mechanism in the form of a ML procedure to do the regression: this would be akin to viewing these as data that are both left- and right-censored. – whuber Feb 06 '13 at 21:36
  • Try beta distribution with parameters smaller than unity, http://en.wikipedia.org/wiki/File:Beta_distribution_pdf.svg – Alecos Papadopoulos Jan 09 '14 at 21:26
  • 3
    This type of bathtub or u-shaped distribution is common in magazine readership where many people will read a single issue of a publication, e.g., in a doctor's office or else are subscribers who see every issue with a smattering of readers in-between. Several comments and responses have pointed to the beta distribution as one possible solution. The literature I'm familiar with points to the beta-binomial as the better fitting option. – user78229 Sep 26 '16 at 17:56

3 Answers3

49

Methods of censored regression can handle data like this. They assume the residuals behave as in ordinary linear regression but have been modified so that

  1. (Left censoring): all values smaller than a low threshold, which is independent of the data, (but can vary from one case to the other) have not been quantified; and/or

  2. (Right censoring): all values larger than than a high threshold, which is independent of the data (but can vary from one case to the other) have not been quantified.

"Not quantified" means we know whether or not a value falls below (or above) its threshold, but that's all.

The fitting methods typically use maximum likelihood. When the model for the response $Y$ corresponding to a vector $X$ is in the form

$$Y \sim X \beta + \varepsilon$$

with iid $\varepsilon$ having a common distribution $F_\sigma$ with PDF $f_\sigma$ (where $\sigma$ are unknown "nuisance parameters"), then--in the absence of censoring--the log likelihood of observations $(x_i, y_i)$ is

$$\Lambda = \sum_{i=1}^n \log f_\sigma(y_i - x_i\beta).$$

With censoring present we may divide the cases into three (possibly empty) classes: for indexes $i=1$ to $n_1$, the $y_i$ contain the lower threshold values and represent left censored data; for indexes $i=n_1+1$ to $n_2$, the $y_i$ are quantified; and for the remaining indexes, the $y_i$ contain the upper threshold values and represent right censored data. The log likelihood is obtained in the same way as before: it is the log of the product of the probabilities.

$$\Lambda = \sum_{i=1}^{n_1} \log F_\sigma(y_i - x_i\beta) + \sum_{i=n_1+1}^{n_2} \log f_\sigma(y_i - x_i\beta) + \sum_{i=n_2+1}^n \log (1 - F_\sigma(y_i - x_i\beta)).$$

This is maximized numerically as a function of $(\beta, \sigma)$.

In my experience, such methods can work well when less than half the data are censored; otherwise, the results can be unstable.

Here is a simple R example using the censReg package to illustrate how OLS and censored results can differ (a lot) even with plenty of data. It qualitatively reproduces the data in the question.

library("censReg")
set.seed(17)
n.data <- 2960
coeff  <- c(-0.001, 0.005)
sigma  <- 0.005
x      <- rnorm(n.data, 0.5)
y      <- as.vector(coeff %*% rbind(rep(1, n.data), x) + rnorm(n.data, 0, sigma))
y.cen           <- y
y.cen[y < 0]    <- 0
y.cen[y > 0.01] <- 0.01
data = data.frame(list(x, y.cen))

The key things to notice are the parameters: the true slope is $0.005$, the true intercept is $-0.001$, and the true error SD is $0.005$.

Let's use both lm and censReg to fit a line:

fit <- censReg(y.cen ~ x, data=data, left=0.0, right=0.01)
summary(fit)

The results of this censored regression, given by print(fit), are

(Intercept)           x       sigma 
  -0.001028    0.004935    0.004856

Those are remarkably close to the correct values of $-0.001$, $0.005$, and $0.005$, respectively.

fit.OLS <- lm(y.cen ~ x, data=data)
summary(fit.OLS)

The OLS fit, given by print(fit.OLS), is

(Intercept)            x  
   0.001996     0.002345

Not even remotely close! The estimated standard error reported by summary is $0.002864$, less than half the true value. These kinds of biases are typical of regressions with lots of censored data.

For comparison, let's limit the regression to the quantified data:

fit.part <- lm(y[0 <= y & y <= 0.01] ~ x[0 <= y & y <= 0.01])
summary(fit.part)

(Intercept)  x[0 <= y & y <= 0.01]  
   0.003240               0.001461

Even worse!

A few pictures summarize the situation.

lineplot <- function() {
  abline(coef(fit)[1:2], col="Red", lwd=2)
  abline(coef(fit.OLS), col="Blue", lty=2, lwd=2)
  abline(coef(fit.part), col=rgb(.2, .6, .2), lty=3, lwd=2)
}
par(mfrow=c(1,4))
plot(x,y, pch=19, cex=0.5, col="Gray", main="Hypothetical Data")
lineplot()
plot(x,y.cen, pch=19, cex=0.5, col="Gray", main="Censored Data")
lineplot()
hist(y.cen, breaks=50, main="Censored Data")
hist(y[0 <= y & y <= 0.01], breaks=50, main="Quantified Data")

Plots

The difference between the "hypothetical data" and "censored data" plots is that all y-values below $0$ or above $0.01$ in the former have been moved to their respective thresholds to produce the latter plot. As a result, you can see the censored data all lined up along the bottom and top.

Solid red lines are the censored fits, dashed blue lines the OLS fits, both of them based on the censored data only. The dashed green lines are the fits to the quantified data only. It is clear which is best: the blue and green lines are noticeably poor and only the red (for the censored regression fit) looks about right. The histograms at the right confirm that the $Y$ values of this synthetic dataset are indeed qualitatively like those of the question (mean = $0.0032$, SD = $0.0037$). The rightmost histogram shows the center (quantified) part of the histogram in detail.

gung - Reinstate Monica
  • 145,122
whuber
  • 322,774
  • great answer (+1). If we were to visually remove the two censoring spkies, it seems to me that the dependent variable has something close to an exponential distribution, as if the underlying data were the length of some process. Is this something to take unto account? – user603 Feb 07 '13 at 00:01
  • @user603 I simulated the quantified values with part of the upper arm of a Gaussian, actually :-). We have to take care here, because the relevant probability model pertains to the residuals and not to the response variable itself. Although it's a little tricky, one can make censored residual plots and even censored probability plots to assess goodness of fit to some hypothetical distribution. – whuber Feb 07 '13 at 00:04
  • my point is that with a doubly censored Gaussian, the histogram of he uncensored values should be somewhat flat, but they seem to be declining gently as we move away from 0. – user603 Feb 07 '13 at 00:10
  • 1
    @user603 Ah, no, that's not the case: take a look at the histogram of the quantified values yourself. They will appear to slope almost linearly downward, exactly as in the question. – whuber Feb 07 '13 at 00:13
  • 3
    I tried censored regression on my dataset and the results crossvalidated better than those from OLS did. A nice addition to my toolkit - thx. – rolando2 Feb 08 '13 at 17:08
10

Are the values always between 0 and 1?

If so you might consider a beta distribution and beta regression.

But make sure to think through the process that leads to your data. You could also do a 0 and 1 inflated model (0 inflated models are common, you would probably need to extend to 1 inflated by your self). The big difference is if those spikes represent large numbers of exact 0's and 1's or just values close to 0 and 1.

It may be best to consult with a local statistician (whith a non-disclosure agreement so that you can discuss the details of where the data come from) to work out the best approach.

Greg Snow
  • 51,722
  • 2
    "Thinking through the process" is good advice. Although beta is a tempting model based on a qualitative inspection of the histogram, I think that if you look closely at the values strictly between $0$ and $1$ you will find they depart substantially from any beta distribution. – whuber Feb 06 '13 at 21:38
  • Yes, always in the range of 0% to 1%...And these spikes are indeed exactly at 0% and 1%. Are zero- or zero-and-one-inflated models applicable to non-count data like these? – rolando2 Feb 06 '13 at 21:43
  • 1
    There is a zero inflated normal, but that doesn't apply here. – Peter Flom Feb 06 '13 at 21:45
  • With exact 0's and 1's but continous between it looks like there may be an underlying continuous distribution with vaues outside that region rounded to 0 or 1. This would be a doubly censored case and models could be fit using that idea. – Greg Snow Feb 11 '13 at 21:14
5

In concordance with Greg Snow's advice I've heard beta models are useful in such situations as well (see a Smithson & verkuilen, 2006, A Better Lemon Squeezer), as well as quantile regression (Bottai et al., 2010), but these seem like so pronounced floor and ceiling effects they may be inappropriate (especially the beta regression).

Another alternative would be to considered types of censored regression models, in particular the Tobit Model, where we consider the observed outcomes to be generated by some underlying latent variable that is continuous (and presumably normal). I'm not going to say this underlying continuous model is reasonable given your histogram, but you can find some support for it as you see the distribution (ignoring the floor) has a higher density at lower values of the instrument and slowly curtails to higher values.

Good luck though, that censoring is so dramatic it is hard to imagine recovering much useful information within the extreme buckets. It looks to me like nearly half of your sample falls within the floor and ceiling bins.

Community
  • 1
Andy W
  • 16,026