In my project, we performed a RNAseq of pig embryos, whose mothers were supplemented or not with arginine on diet. We are trying to analyse the differential expression of this data using DESeq2 with a nested model: Y= mean + Sex + treatment + Sex by treatment interaction + mother within treatment + error , but it is not working, because mother should be considered as a random effect.

When we tried this model with fixed effect, the error message "Error in checkFullRank(full) : the model matrix is not full rank, so the model cannot be fit as specified. Levels or combinations of levels without any samples have resulted in column(s) of zeros in the model matrix. Please read the vignette section 'Model matrix not full rank. Vignette('DESeq2')" it is returned to us.
We have followed the suggestion but we still not have success.

We have been extensively looking the DESeq2 manual, papers and forums to solve this issue, but we haven't found any similar example available. Therefore, we would like to know if it is possible to use the DESeq2 analysis with random effects?

Hi, everyone.
I have performed my in Limma package as suggested because there are nested effect interaction (random effect) in statistical model.
Looking to my results the log foldchange seems to be wrong to me. Some genes have too high LogFC values such as -14997.66667 with p-adjust 0.02518022, for example.
We don't have experience with this package, we would like to know if this is normal.

To my knowledge, none of the negative binomial GLM-based packages (e.g. edgeR & DESeq2) support models with random effects. Limma is the only package I know of that does support RNA-seq analysis with random effects. It supports a model with any number of fixed effects and one factor as a random effect, which seems to match your use case. The Limma User's Guide shows how to do this using the duplicateCorrelation function.

SAMPLE_ID Mae Tissue Sexo Treatment
P8175-25DA-E1M P-8175 embriao M ARG
P8175-25DA-E2F P-8175 embriao F ARG
P8175-25DA-E5F P-8175 embriao F ARG
P8191-25DA-E1F P-8191 embriao F ARG
P8191-25DA-E2F P-8191 embriao F ARG
P8191-25DA-E5M P-8191 embriao M ARG
P8192-25DC-E1M P-8192 embriao M CONT
P8192-25DC-E2F P-8192 embriao F CONT
P8192-25DC-E5M P-8192 embriao M CONT
P8193-25DA-E1F P-8193 embriao F ARG
P8193-25DA-E2F P-8193 embriao F ARG
P8193-25DA-E5F P-8193 embriao F ARG
P8206-25DC-E1M P-8206 embriao M CONT
P8206-25DC-E2F P-8206 embriao F CONT
P8206-25DC-E5F P-8206 embriao F CONT
P8285-25DC-E1F P-8285 embriao F CONT
P8285-25DC-E2M P-8285 embriao M CONT
P8285-25DC-E5M P-8285 embriao M CONT

Yes, essentially you have the same setup as the example above. You want to control for the mother, which is nested within group (here, treatment). You want to estimate the sex effect, controlling for mother, and see if it is different across treatment groups. You can read the vignette section and substitute treatment for "group" and sex for "condition". The one important distinction is that in DESeq2, you can't estimate the treatment effect and simultaneously control for mother, because those are directly confounded. You would have to use the duplicateCorrelation approach to also perform that comparison.

Hi, everyone. I have performed my in Limma package as suggested because there are nested effect interaction (random effect) in statistical model. Looking to my results the log foldchange seems to be wrong to me. Some genes have too high LogFC values such as -14997.66667 with p-adjust 0.02518022, for example. We don't have experience with this package, we would like to know if this is normal.

I moved your post from an "Answer" to a "Comment", because it doesn't answer your original question.

Also, I'd recommend that you formulate a new post, tagged with "limma", if you want to ask specific questions of the package authors.