I wanted to ask if using the following model design was appropriate for the question I'm trying to answer. I have samples that are paired (BG and AG status). The samples come from 8 individuals and each individual has a BG and AG sample. Additionally, there are 4 treatment groups (A, C, T, E) and a control (W).
I wanted to compare the difference between the two statuses across the treatments and control. I was able to use this model in DESeq2, but my local statistician suggested that I compare the results to a package that uses a different distribution.
DESeq2 model: ~treatment + treatment:status
Based on this site, I created a similar model, but received this output:
# Define the metadata categories TB.treatment <- pData(n_paired_merged_MRE)$Treatment TB.ID <- pData(n_paired_merged_MRE)$id TB.bs <- pData(n_paired_merged_MRE)$BrushStatus # Define the normalisation factor normfactor_MRE = normFactors(n_paired_merged_MRE) normfactor_MRE = log2(normfactor_MRE/median(normfactor_MRE)+1) # Create the model TB.mod <- model.matrix(~ TB.treatment + TB.bs + TB.ID + normfactor_MRE) settings <- zigControl(maxit = 10, verbose = TRUE) TB.fit <- fitZig(obj = n_paired_merged_MRE, mod = TB.mod, useCSSoffset = FALSE, control = settings) # output Coefficients not estimable: TB.IDW21AG TB.IDW21BG normfactor_MRE TB.IDC21BG TB.IDE21BG TB.IDT21BG it= 0, nll=NaN, log10(eps+1)=Inf, stillActive=1141
Coefficients not estimable: TB.IDW21AG TB.IDW21BG normfactor_MRE TB.IDC21BG TB.IDE21BG TB.IDT21BG Error in lm.fit(mmZero, qlogis(pi)) : NA/NaN/Inf in 'y' In addition: Warning messages: 1: Partial NA coefficients for 1141 probe(s) 2: Partial NA coefficients for 1137 probe(s)
Thanks for the help!