Family effects in DESeq2/limma/nlme

[Lohman]

Hello All,

I'm looking for your opinions on proper analysis of gene expression data from RNAseq with a rather complicated design. I am running a controlled lab infection experiment, exposing two populations of hosts to parasites, and looking for changes in gene expression associated with either infection, or failure of infection (Exposed). Here's the basics:

Total n = 95

Population:
GG:55
RR:40

Family:
19 total
~5 individuals/family (family member are either control, exposed, or infected)

Status:
Control:36
Exposed:37
Infected:22

Batch is the lane they were sequenced on (yes, I should have multiplexed properly, but that mistake has already been made).

I also have data on sex, size of host, and size of parasite.

In the PCA, families cluster together, so I believe that family should be considered a random effect. After looking around a bit on various threads, it's often suggested to use a fixed effect instead of a random effect, but there are usually many fewer levels that I have here (3-4 vs 19).

How would you incorporate family into the analysis? I'm open to all suggestions.

Alternatively, a general strategy for model testing in gene expression would be equally welcome.

Thanks in advance for you help

Lohman

[dpryan]

The number of levels won't matter, just use a fixed effect with something like DESeq2 that's already doing shrinkage anyway. This is identical to every clinical RNAseq experiment published, where there's a level for every patient.

[Lohman]

Hi Devon,

Thanks for your reply. If you fit family as a fixed effect in DESeq2 then you must drop population from the model because family predicts population. The resulting model would identify differences in infection status (treatment) but ignore differences in population. Would you attempt to interpret the family predictor at all? Or simply ignore it?

Thanks for your suggestions.

Lohman

[dpryan]

If all members of a family are in a single population then population becomes the average of the family effect. Normally one doesn't care about family effects, since it's just something that needs to be compensated for.

[Lohman]

Would you then run two different models, one with population and the other with family and use them for two separate purposes (including interactions in the model with population)?

[dpryan]

If you actually have a notable family effect (not unusual) and your population is stratified/partitioned by family then you won't get a very interpretable population effect regardless of the model or method you use. If that's not the case, then just use a contrast (yes, this can get annoying very quickly, though you might be able to just include the interaction in the original model).