Hi all,

I have the following data:

50 subjects that belong to 2 groups (25 in each group). For each subject I
have replicated measurements (in the range of 15-30 replicates per subject)
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 in
each replicate is assumed to be a binomially distributed random variable.
To complicate things even more, for each subject from each group I have not
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 the
success probabilities of subjects from group G2.

In other words, for each subject in each replicate a success probability
can be estimated from the observed data. Then the variance of those success
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 itself
   affects the accuracy of estimating the success probability, but in
   addition, I cannot rule out an inherent dependency between the number of
   observed trials and the success probability (there is no control over the
   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.


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 all
replicates, as well as the mean number of trials over all replicates (added
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 + (1|environment),
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.

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 since
   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 normal
   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 = n-1)
   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]]

