[R-sig-ME] A mixed effects model for variances of binomialprobability as a responsee
rdiaz02 at gmail.com
Fri Mar 28 03:06:52 CET 2014
On Thu, 01-01-1970, at 01:00, nimrod.rubinstein <nimrod.rubinstein at gmail.com> wrote:
> Hi Ramon (and others reading this),
> Thanks for the response.
> I'm attaching my real data - it has these columns:
> subject - subject name (S1, S2,...)
> group - group name (G1, G2) - a subject can only belong to either of the
> two groups
> environment - environment name (E1, E2, E3).
> the next 27 columns (named n.1-n.27) are the number of trails in each
> the following 27 columns (named k.1-k.27) are the number of successes out
> of the number of trials in each corresponding replicate (so k.1 is the
> number of success out of n.1 trials, and so on).
> Note that there are NA values in corresponding n and k cells in the table.
> Also, not all subjects were sampled from all 3 environments due to
> technical problems, so the data are not fully balanced.
> I know this table is a bit different from my earlier description of the
> data. It is so only and it is because I wanted to make it simple.
> Let me try and clarify, and hopefully all of the below will answer your
> What may affect the variances of the success probabilities (proportions)?
> 1. The proportion itself (natural to binomially distributed data). See
> attached Fig1 that plots sample variances of success probabilities vs.
> their sample mean.
But what are those variances? How are you computing them? (See below, though)
> 2. The number of trails in each replicate. As I said, the number of trails
> is inherent to the subject and the environment it was measured in, the
> experimenter has no control over it. As the attached Fig2 shows, subjects
> with a high average number of trials have lower sample variances of their
> success probabilities. The opposite however is not clear. The sample
> variances of subjects with a low average number of trials vary greatly. So
> obviously, the number of trails has an effect on the variance of the
> proportion but there's still a sizable fraction of variance (of the
> proportion variances) that are not explained by that.
> 3. The group - that's the effect I'm interested in. Fig1 and Fig2 show that
> the G1 dots have lower sample variances of success probabilities than G2
> 4. The environment. These environments were chosen at random and therefore
> I think they should be REs. (In all attached figures the environments isn't
> I've also attached a figure of the success probabilities vs, the number of
> trials (Fig3).
Yes, I think the scenario I was concerned about is not here. Good.
> I think the figures show pretty clearly that both number of trails, success
> probability, and group affect the variance of success probability, and
> therefore the model I'm looking for should estimate all of these.
It is way too late here, but maybe we are using the same name (variance)
for two different things here: the variance of p, which is just a simple
function of the proportion, of the p itself, and the variance in the
collection of p's, where each p is not necessarily estimating the same
probability. And I think you care about the second problem, something like
saying that individuals in G1 and G2 have different dispersion in their
p's. So you have a distribution for the p's of G1 and a distribution for
the p's of G2, and you want to know whether those distributions have the
same variance while also possibly controlling for other factors. I think
the book by Gelman and Hill has some examples of this kind (I do not have
it here now, though).
> Thanks a lot,
> On Thu, Mar 27, 2014 at 4:39 AM, Ramon Diaz-Uriarte <rdiaz02 at gmail.com>wrote:
>> Dear Nimrod,
>> Just my 2 cents (even less than 2):
>> On Thu, 01-01-1970, at 01:00, nimrod.rubinstein <
>> nimrod.rubinstein at gmail.com> wrote:
>> > Hi all,
>> > I have the following data:
>> > 50 subjects that belong to 2 groups (25 in each group). For each subject
>> > have replicated measurements (in the range of 15-30 replicates per
>> > of the number of trails as well as the number successes out of these
>> > trails. Each single trial in each single replicate for each subject is
>> > assumed to be a Bernoulli experiment, therefore the number of successes
>> > each replicate is assumed to be a binomially distributed random variable.
>> > To complicate things even more, for each subject from each group I have
>> > 1 but 3 series of these replicated measurements, from 3 different
>> > environments (the 3 environments are the same for all subjects).
>> > A few example lines of how my data looks like:
>> > subject group environment number.trails number.successes
>> > S1 G1 E1 10100,13209,17621 8500,11622,14684
>> > S1 G1 E2 13322,12225,10024 10351,10222,7951
>> > S1 G1 E3 9812,13500,13214 7245,10622,10684
>> > .
>> > .
>> > .
>> > S50 G2 E1 60547,75345,68412 12014,13214,15244
>> > S50 G2 E2 50897,62541,61247 10141,13251,30645
>> > S50 G2 E3 59214,84264,74156 17245,16425,11123
>> > Where the comma separated number of successes correspond to the comma
>> > separated number of trials in each replicate.
>> > I want to test whether the variances of the success probabilities of
>> > subjects from group G1 are significantly smaller than the variances of
>> > success probabilities of subjects from group G2.
>> But, as you say in point 1. below (and in another reply), if this is really
>> binomial, then you already know that if probabilities in G1 are different
>> from probabilities in G2 your variances will differ.
>> So I am not sure I even understand this hypothesis/question.
>> > In other words, for each subject in each replicate a success probability
>> > can be estimated from the observed data. Then the variance of those
>> > probabilities over all replicates from a specific environment can be
>> > estimated. I would like to know whether these variances are smaller in
>> > subjects from G1 vs. subjects from G2.
>> > The issues are:
>> > 1. The variance of a success probability depends on the success
>> > probability itself (natural to binomially distributed data). I need to
>> > account for this effect.
>> > 2. The number of trails (although at least in the hundreds for all
>> > subjects) varies several orders of magnitude among subjects. That by
>> > affects the accuracy of estimating the success probability, but in
>> > addition, I cannot rule out an inherent dependency between the number
>> > observed trials and the success probability (there is no control over
>> > number of trials in the experiment - it is observed). In other words,
>> > having a lower number of trails by itself may affect the success
>> > probability.
>> That last issue seems very problematic to me. If I understand correctly, we
>> could think of a perverse situation where, say, you can never observe low
>> success probabilities if those are related to number of trials around 0.
>> Can you do a simple plot of probability of success (or its logit) vs number
>> of trials, to rule out perverse scenarios? Is probability of success always
>> larger than 0? And number of trials always sufficiently different from 0?
>> And you say "varies by subject", so the number of trials is always the same
>> for a subject? And the potential relationship prob <-> number of trials is
>> due to among-subject variation? So you could approach the problem via a
>> random effect for subject? If you compare G1 and G2 for number of trials:
>> are there differences?
>> > I though a simple way to test my hypothesis is to compute the success
>> > probability in each replicate as number.successes/number.trails, and then
>> > the sample mean and sample variance of these success probabilities over
>> > replicates, as well as the mean number of trials over all replicates
>> > to my data as sample.mean, sample.var and mean.trials fields). Then I can
>> > use the lmer function in the lme4 package this way:
>> > null.model = lmer(sample.var ~ sample.mean + mean.trials +
>> > my.data, REML = F)
>> > alt.model = lmer(sample.var ~ sample.mean + mean.trials + group +
>> > (1|environment), my.data, REML = F)
>> > So fixed effects are: sample.mean, mean.trials and group, and random
>> > effects is environment since the three environments in my data are
>> simply a
>> > random sample from many other possible environments.
>> I do not see how this would answer your questions, given the above
>> concerns/issues and other possible approaches (just compare prob. between
>> G1 and G2, accounting for environment).
>> > To test whether group has a significant effect I'll do:
>> > anova(null.model, alt.model).
>> > My concerns are:
>> > 1. Should I add subject as a random effect? My intuition is not to
>> > subjects are unique to groups, which is the effect I'm after.
>> However, I do
>> > want to account for the fact that sample variances of the same
>> subject from
>> > the three different environments are correlated.
>> > 2. The response - sample variances, cannot be assumed to follow a
>> > distribution (they only take positive values to begin with and as far
>> as I
>> > know (n-1)*sample.variance/population.variance ~ chi.square with df =
>> > and therefore the assumptions of a linear mixed effects model do not
>> > strictly hold. Is there a generalized linear model that would be
>> > appropriate for this?
>> > Thanks a lot,
>> > nimrod
>> > [[alternative HTML version deleted]]
>> > _______________________________________________
>> > R-sig-mixed-models at r-project.org mailing list
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>> Ramon Diaz-Uriarte
>> Department of Biochemistry, Lab B-25
>> Facultad de Medicina
>> Universidad Autónoma de Madrid
>> Arzobispo Morcillo, 4
>> 28029 Madrid
>> Phone: +34-91-497-2412
>> Email: rdiaz02 at gmail.com
>> ramon.diaz at iib.uam.es
Department of Biochemistry, Lab B-25
Facultad de Medicina
Universidad Autónoma de Madrid
Arzobispo Morcillo, 4
Email: rdiaz02 at gmail.com
ramon.diaz at iib.uam.es
More information about the R-sig-mixed-models