I am in search of a fixed intercept (through the origin) random slopes model. However, instead of additive errors, I want the SLOPE to be the source of random error.

The problem is this: for a known volume of fluid, we are performing a regression between the volume and the surface area created when that volume of fluid is dropped onto a surface. Clearly, when the volume is zero, the area will be zero, and measurements were taken very close to zero (motivating regression-through-the-origin). It is also reasonable to assume that each unit has a different mean slope (unique to that unit), so we can think of each unit's average slope as a sample taken from a population with an overall mean slope. Furthermore, I want to show that the mean slope for the populations increases when the fluid in diluted.

10 units were randomly selected. Each had fluid removed, then that fluid was separated and one partition was diluted. A fixed volume was drawn from each partition and dropped onto an experimental surface (assume no surface effect). This process was repeated 12 times for each volume in each partition in each patient. The surface areas were then measured.

To summarize the problem:

For Unit $i$, $i \in 1,...,10$

Treatment $j$, $j \in \{0,1\}$. $0 = $ undiluted, $1 = $ diluted

Volume $k$, $k\in 1,...,7$,

$v_k \in \mathbb{R}^{+}$, taking fixed values $ \{ 2.5, 5, 7.5, 10, 12.5, 15, 20 \}$

Replicate $l$, $l \in 1,...,12$


  • $\beta_{ij}$ ~ $N(\mu_{j},\sigma^2_{\beta})$ , i.e. $\beta_{ij}=\mu_{j} + \gamma_i$. Each randomly selected unit has it's own characteristic slope, from a population whose mean slope depends on the treatment.
  • $B_{ijk}$ ~ $N( \beta_{ij} , \sigma^2_e )$, i.e. $B_{ik}=\beta_i + \varepsilon_{k}$
  • $A_{ijk} = B_{ik}v_j$

I tried simulating a dataset and running lme4 on it. My code to simulate is:

rm(list = ls())
I <- 10 # 10 Units
J <-  2 #  2 treatment groups
K <-  7 #  7 discs, one for each volume
  v = c(2.5, 5, 7.5, 10, 12.5, 15, 20)
L <- 12 # 12 replicates for each volume on a single disc

test.df <- data.frame( Units = as.factor(sort(rep(c(1:I),J*K*L))),
                      Trt = rep( sort( rep( 0:1 , K*L )) , I ), 
                      Disc = rep(sort(rep(c(1:K),L)), I*J), 
                      Volume = rep( sort(rep( v , L )) , I*J ),
                      Rep = rep( 1:L , I*J*K ) );head(test.df,18)

# the slope is a N(14,0.5) for undiluted.  N(16,0.5) for diluted
beta <- 14 + 0.5*rnorm(I)
trt.increase.beta = 2*test.df$Trt
test.df$beta <- beta[test.df$Patient] + trt.increase.beta

# The observed data comes from (fixed volume)*(the unit's slope + zero mean replication error with variance 0.25)
test.df$Area <- test.df$Volume * ( test.df$beta + 0.5*rnorm(1) )

plot(test.df$Volume , test.df$Area, pch = 20)

How to I model this in R??? What is the syntax for the lme4 package? I have searched in several papers and blogs, but haven't come across an examples.

What I have now is:

random_slope_through_origin.lm <- lmer( Area ~ 0 + Volume + Trt*Volume + ( 0 + Rep | Volume ) + ( 0 + Units | Volume ) , data = test.df)

When I run this, I get all kinds of errors, most of which I believe comes from the heteroscedasticity of the data (the variance increases with the volume). I think lme4 is modeling additive variance. Is there a way to measure multiplicative variance?

Thanks for any help you might offer.

  • $\begingroup$ Surface areas are necessarily positive (or at the very least, non-negative, if you have actual zero-volumes); it doesn't make sense to me to assume they're normal. Their logs might be closer to normal. Or you might take some other positive observations with multiplicative noise type of model. $\endgroup$
    – Glen_b
    Jun 28, 2017 at 6:38
  • $\begingroup$ True, but you see I am assuming that the SLOPE is normal. A normal random variable multiplied by zero will be zero. Furhtermore, a normal random variable with a positive mean and small variance is very unlikely to ever become negative. Also, while I agree that some other non-negative distribution may be more realistic, this is true or MANY datasets, yet it is still common practice to assume approximate normality and proceed accordingly. $\endgroup$
    – Lewkrr
    Jun 28, 2017 at 16:26
  • $\begingroup$ Is normality necessary for this? I was trying to point you to potentially simpler analyses that would be easier to justify. Sure, normality might be a reasonable approximation, but there may be other good approximations that work at least as well and may be fairly simple to deal with. If you want normality, no problem, go right ahead, but sometimes people don't consider anything else, and it's not always as good as some other choice. $\endgroup$
    – Glen_b
    Jun 28, 2017 at 23:40
  • $\begingroup$ @Glen, thanks for the feedback! It's greatly appreciated.... I have considered a gamma regression, with the shape = $\frac{\beta^2}{\gamma}x^{1-\tau}$ and scale = $\frac{\gamma}{\beta}x^{\tau}$ so the $E(Area) = \beta x$ and the $Var(Area) = \gamma x^{1+\tau}$, $\tau \in \left[ -1 , \infty \right)$, however I encounter convergence issues. What kind of simpler analysis do you propose? How would you approach this? $\endgroup$
    – Lewkrr
    Jun 29, 2017 at 19:15


Your Answer

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