In this question/thread I'd like to hear your opinions about various designs of typical cell culture experiments. It's a complex matter and I have several questions but I'll try to squeeze most of it in this opening post. For model designs I'll use the lme4 notation, with (1|factor) denoting random factors, I'm rather new to this, so feel free to point out errors in interpretation you feel I might be making. I've also made a sort of example sheet to illustrate some of my points, as I'm not versed well enough in the subject to formulate text unambiguously.
The type of experiments I mean are your run-of-the mill micro-plate experiments, with a normally distributed continuous dependent variable, let's call it DV (think for instance absorption on a plate reader). A typical design might look as follows: a (96-well) microplate is seeded with cells and the cells are subjected to a few experimental conditions. Several wells are exposed to the same condition (these are called technical replicates here), the number of which per condition usually depends on how many wells are left on the plate and aesthetic preferences in plate layout design. Note: adding extra wells per condition is very easy/non-time consuming.
Traditionally the technical replicates are averaged, which makes the plate (or perhaps batch if there are several plates) the experimental unit and the experiment is repeated on different days. The formal reason to use averages is that values within a plate *might* be (auto-)correlated and it is difficult to account for this using traditional/frequentist methods.
It is however possible to correctly model the full design using mixed models, by for instance analysing the data as nested within plate. So now to my questions:
A: Do you think it is valid in general to use mixed models in this situation?
A1: Papers with in vitro studies have poor reproducibility, will using more powerful models further increase Type 1 errors?
B: Would lmer(DV~Condition + (1|plate)) be the correct way to specify the model? I guess plate is a crossed factor here, do you have to make nesting more explicit somehow? (see also sheet)
C: Does the number of technical (= lowest level) replicates influence the power of the experiment much? Note again that adding extra technical replicates in the design is very easy and can go up to large numbers. If it does, how does this affect your type I error? BTW: I realize increasing technical replicates will increase precision of the estimate, even if you do reduce them to an aggregate.
C1: How does the actual amount of autocorrelation affect the power? Alternatively put: how does the ratio of the within-plate variance and the between-plate variance affect power?
C2: If you said 'no' at A, doesn't the relationship between the variance at replicate and at plate level feel like important information somehow? Shouldn't it be in your ideal model?
D: Is there any difference in interpretation of the results compared to the traditional aggregated model (i.e. lm (aggrDV~Condition))
D1: Will the mixed model always be more powerful than the model with aggregated data? Does it approach power of the aggregated model when lowest level variance is comparatively low or comparatively high?
Now a special case that pertains to my own experiments: I'm doing experiments with primary cells from diseased donors, which are rare to come by. In this case I would traditionally use 'donor' as experimental unit (note in the first part all cells are clones from the same tumor/donor). As donors are rare I'm having power issues in my design (low n, no way of increasing it).
E: If you said 'no' at A: Do you think a design where you repeat experiments with the same donor several times and nest plate within donor is valid?
E1: Is lmer(DV~Condition + (1|donor) + (1|plate) then correct? Does order matter? Would you take 'measured values' or 'DV' as dependent variable? (see 2nd sheet)
E2: If your said no at E, won't that mean that you'll have to discard almost every paper using a cell-line?
F: In the simulation sheet donor example I've assumed plate is nested within donor, but in reality you just have the measured values. If you don't do a hierarchical model it starts to matter whether you aggregate by donor or by plate. How do you know which is correct? (see Aggr. DV1 vs Aggr DV2).
A: yep, I think that this is the way to go -- at least I do use mixed models in this situation.
A1: I am a biologist by training, so my natural instinct is to plan an alternative experiment or increase sample size rather than fiddle too much with model intricacies. However in the particular set up you are describing, I think that averaging the values makes you loose information. Before rushing into the complexity of different models, I would nonetheless first try to see where the main sources of variability come from.
For example, I had to analyse data from an experiment in which three biological replicates (independently differentiated cells) were analysed in three experiments, each with three replicates (independent infections, let's call them "technical" although they were not really only technical replicates). While there was some variability within an experiment, it was nothing compared to the effect of the biological replicates, and in the end both averaging out the technical replicates or a mixed model yielded almost identical results.
Formally, you would create different models, with and without the within-experiment variability, and compare them with ANOVA. See the Pinheiro/Bates book for details.
B: looks good to me. I use the nlme package in R, though (and follow the Pinheiro book in most details).
C: I always argue for increasing the biological replicates rather than the technical replicates. You see, increasing the number of technical replicates improves the answer to question "what is the response in this particular biological replicate", while usually you would like to answer the question "what is the response in general, common for biological replicates".
C1: no idea.
E1 (and other questions as well): I would try different designs and check whether the models are significantly different, and if no, go for the simplest design.
Thanks Marianne and January! Glad somebody took the time to read all of my questions, I think I might have gone a bit overboard with the number of questions.
Marianne: You're probably right that GEEs are a simpler solution in many cases. To be honest I don't really understand the difference between GEEs and mixed models (to me they are similar in application). From what I understand (in my layman's terms) GEEs ignore the covariance structure in the model somehow and basically set your random factors to 1. I'm not sure if the end results differs much with continous data, perhaps the confidence intervals are more robust? (This stuff is out of my league).
In the articles citing the paper you supplied, I did find some articles that use mixed models which actively take technical replicates into account. For example this one: http://linkinghub.elsevier.com/retrieve/pii/S0003269712002576
In this article they do find smaller S.E.s in the mixed model compared to the aggregated/averaged model, thus showing the approach has a bit more power. This artice: http://linkinghub.elsevier.com/retrieve/pii/S0888754309000986 actually looks at type I error rate and finds it to be at the nominal alpha level.
BTW: I did consult a statistician for the design and the verdict was that mixed models are (mathematically at least) appropriate in this case. As for the interpretation and which conclusions could be draw, she couldn't comment however.
January: Good to know there's some more experience with mixed models out there. Indeed the variance of technical replicates is usually minimal (which is exactly why it feels like cheating/inflating n when you include them in the model, even though you 'waste data'). If it doesn't matter (much) in the final outcome, I paradoxcially feel safer to include them. Besides. the way you aggregate (mean/median) can matter quite a bit in some cases, and with a hierarchical model you don't run the risk of making an interpretation error at this point.
I'll try and run some simulations for C and C1 if it isn't immediatelly apparent, it might help in the development of some intuition in model design (although increasing biological replicates, as you say, is inherently better just from the conceptual point of which variance you are addressing).
As for E, I'm stuck with my limited donor numbers, but I think it's valid in either case to measure each donor (derived celline) repeatedly to boost my power a bit even though it adds a lot of complexity to the design.
So, I played around a bit, trying to get the hang of this stuff.
I failed :-S
Attached is a file with some lines of code where I tried to make some example data from scratch and then compare different ways of analysis. The data looks plausible enough and I get expected Type I error (like-values) with the t-test. When I use mixed models however, I get very weird results. There's probably something very stupid going on, but I'm not seeing it. Help?
P.S. If you find yourself facepalming at my n00b code, feel free to upload amended versions.
Ugh, without now going into details (your code is no worse than mine, but in reading code from anyone else it is nice to have *plenty* of comments), I see that your plates are nested within conditions; you are calculating the per-plate averages independently for conditions A and B. Which means that you should not analyse the interaction Plate * Condition, since "plate 1" in condition A is not the same as "plate 1" in condition B; the plates in B are independent of plates in A.
Either rename the plates:
testdata$platenew
Hi January, I included both nlme and lme4 models to show that they are the same. In this (balanced) case nlme would indeed be preferable. I am aware of (which is not to say I understand completely) the discussion about p-values. In this case of a well-controlled prospective cell-culture experiment you do want to make a a classical inference/hypothesis testing though, so p-values are appropriate (I think).
I've used the Anova() function from car to derive a p-value, for lme class models it gives the same value as anova(), if you specify the same model with lme4, i.e. lmer class it will return an identical value. I've also done checks with MCMC procedures, and the results are very very similar.
To check my understanding of what's going and if I'm correct.
Plates nested within conditions: isn't that what's going on in my experiment? Say I have a (96-)well plate with cells (plate 1), I stimulate 3 wells with a compound (condition B, a cytostatic), and leave 3 well unstimulated (condition A). I do the same on different days (Plate 2:3). Each time I measure a variable (I don't now, say viability by MTT). Aren't then plates nested in conditions?
As for interpretation of the models, I want to correct mainly for seeding density of the cells (which differs per plate, but is assumed equal per wells within plate). This should be the model intercepts per plate, right? i.e. (1|Plate) I might expect cells to react differently per plate (delta viability, depending perhaps on seeding density), which should be 1+Condition|Plate, but I think this might be a minor effect.
Either way, now matter how I specify (1|Plate), (1+Condition|Plate) or both or even recoding for independent plates like you propose, I get a lot more than my nominal alpha if I set both means to the same value with either nlme or lme4 (but not with the t-test, which is at alpha).
Hmm. I'll have to re-read Bates again I guess...
- Badges
- Science topic