Paired Difference Analysis
0
0
Entering edit mode
nhua • 0
@nhua-13743
Last seen 3.5 years ago
University of Texas at Austin

Hello,

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). 

Sample ID Treatment Status
A01BG A BG
A01AG A AG
C01BG C BG
C01AG C AG

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

results(dds, contrast=list("TreatmentC.StatusBG","TreatmentA.StatusBG"))

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! 

metagenomeseq • 470 views
ADD COMMENT

Login before adding your answer.

Traffic: 280 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6