Question: DESeq2 design and stats question
gravatar for rbronste
9 months ago by
rbronste60 wrote:

So will just list current design first and then get into the question:

rangedCounts <- dba.peakset(Adult_BN_count, bRetrieve=TRUE)

nrows <- 1054182 
ncols <- 12
counts <- matrix(runif(nrows * ncols, 1, 1e4), nrows)

sampleNames<-c("MBV1",    "MBV2",    "MBV3",    "FBV1",    "FBV2",    "FBV3", "MBE7", "MBE8", "MBE9", "FBE1", "FBE2", "FBE3")
sampleSex<-c("male", "male", "male", "female", "female", "female", "male", "male", "male", "female", "female", "female")
sampleTreatment<-c("vehicle", "vehicle", "vehicle", "vehicle", "vehicle", "vehicle", "BB", "BB", "BB", "BB", "BB", "BB")
colData<-data.frame(sampleName=sampleNames, sex=sampleSex, treatment=sampleTreatment)

# Retrieve counts from dba.count without the site interval ranges

counts <- as.matrix(mcols(rangedCounts))

# Construct a SummarizedExperiment

se<-SummarizedExperiment(assays=list(counts=counts),rowRanges=rowRanges, colData=colData)

# DESeq2 analysis

table(sampleSex, sampleTreatment)

dds<-DESeqDataSet(se, design= ~sex + treatment)

dds$group <- factor(paste0(dds$sex, dds$treatment))

design(dds) <- ~0 + group

dds <- dds[ rowSums(counts(dds)) > 1, ]

dds <- DESeq(dds, betaPrior = FALSE)

So my question is about two types of results I am aiming for, group comparisons and pairwise comparisons, results method listed below:

#pairwise comparison

BN_group_male<-results(dds, format = c("GRanges"), contrast=c("group", "maleBB", "malevehicle"))

#group comparison

resFemVeh_BN<-results(dds, lfcThreshold = 1, altHypothesis = "greater", format = c("GRanges"), contrast=c(-1/3, 1, -1/3, -1/3))

My overall goal here is to see:

1. Example: What are the male vs female diff irrespective of treatment etc.

2. Example: What are differential peaks specific to only the Male BB condition and no other.

I am wondering if this is the proper way to get the results for these particular comparisons given the above design? 

ADD COMMENTlink modified 9 months ago by Michael Love20k • written 9 months ago by rbronste60
gravatar for Michael Love
9 months ago by
Michael Love20k
United States
Michael Love20k wrote:

The first is straightforward: a comparison of the BB vs vehicle in males. The second is a comparison of the second group against an average of the others. Use plotCounts to get a sense of the top genes that the results tables are pulling out.

ADD COMMENTlink written 9 months ago by Michael Love20k


My central question here is whether the design makes sense for the two ways in which I am retrieving results? 

Also, in addition, I was noticing that when I lfcThreshold the #1 result type (a comparison of the BB vs vehicle in males) at the same cutoff as #2 I get pretty much no diff peaks (hence no lfc thresholding in the above #1 results formula), and wondering why that might be? Thanks again.

ADD REPLYlink modified 9 months ago • written 9 months ago by rbronste60

"Male vs female irrespective of treatment" you could get by averaging male vs female in each group. So putting 1/2 and -1/2 in the contrast. You don't seem to have anything like that in your code.

One vs the rest you can get with your second contrast.

Getting no diff peaks with cutoff: do some exploration and see if you find the answer. I would look at e.g. an MA plot of the results table.

ADD REPLYlink written 9 months ago by Michael Love20k

Ok great thanks will put in that additional comparison, but in terms of the actual contrast I guess it would be the following, for male vs female only:

contrast=c(1/2, -1/3, -1/2, -1/3)

I guess just a little confused in this format for how to specify the comparison, not as much in the vignette. 

Thanks again!

ADD REPLYlink written 9 months ago by rbronste60

Put 1/2 for coefficients associated with M and -1/2 for the coefficients associated with F. This is a bit of a complex contrast because it isn’t represented by any coefficient in the design. There isn’t a M vs F “irrespective” when you are fitting a group for each combination of sex and treatment. You may want to discuss this in more detail with a statistical collaborator. I mostly provide software feedback here, but have limited time for further statistical instruction.

ADD REPLYlink written 9 months ago by Michael Love20k

Thanks again!

ADD REPLYlink written 9 months ago by rbronste60
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 460 users visited in the last hour