Two parts of small intestine from same individuals, It's technical or biological replicates?

Question

[hjsunbio]

Background We have an RNA-seq project on small intestine. The experimental design were like this: Four physiological states, each with three individuals (biological replicates), each individual with two samples (two parts of small intestine), thus 4 physiological states × three individuals × two parts. The reason we included two parts of same individual’s small intestines is that we want to take consideration of tissue heterogeneity, by pre-analysis of these replicates the difference between two parts is actually not significant.

Question Here is my question, when I applied a generalized linear model for the differential expression analysis on samples from two physiological states (e.g. DESeq2), should I treat samples as independent one thus 6 vs 6, or merge two samples from same individuals thus 3 vs 3. For the later one, can I regard it as technical replicate and merge reads from the two part together, or merge counts calculated from the two parts (e.g. collapseReplicates function in DESeq2). Or, should I include "individuals" as a blocking factor in the differential expression analysis.

There is an intensive discussion about biological replicates and technical replicates, like these Technical Or Biological Replicate Which step is better for merging technical replications? https://www.biostars.org/p/304008/#

but I am still confused about my situation, should I treat the two tissue parts from same individuals as biological or technical replicates.

Answer 1

[Asaf]

From DESeq2 Manual, look for Group-specific condition effects, individuals nested within groups. You can prepare a table with individual number and sample number and then treat each sample as individual within a group (group being individual in your case). The coldata would look like (correct me if I'm wrong):

SampleName PhysiologicalState Individual Sample
Sample1  State1  1  1
Sample2  State1  1  2
Sample3  State2  1  1

Then your model will look like model.matrix(~ PhysiologicalState + PhysiologicalState:Individual, coldata) And then go and remove columns with zeros. Please check if this works, I'm not 100% sure.

[hjsunbio]

Thanks a lot and sorry for the delay response. It is really a big hand for me. Actually, this is a big concern of one reviewer of our work. I have read all the materials you provided, below is my table for construct model matrix. I would like to know whether it is right. Another question is how should I code if I want to know genes differently expressed between two Physiological states. e.g. use DESeq2

  ddsMat_compare <- DESeqDataSetFromMatrix(countData = countdata_compare,
                                        colData = design_compare,
                                        design = ~ PhysiologicalState+PhysiologicalState:individual)
ddsruv_deseq = DESeq(ddsMat_compare)
???

or edgeR

y <- DGEList(counts=countdata_compare, group=conditions)
y <- calcNormFactors(y,method="TMM")
y <- estimateDisp(y, design_edge_used, robust=TRUE)
fit <- glmFit(y, design_edge_used,robust=TRUE)
glmLRT(fit, ???)

samples individual PhysiologicalState

4H_S1_P1 1 4H

4H_S1_P2 1 4H

4H_S2_P1 2 4H

4H_S2_P2 2 4H

4H_S3_P1 3 4H

4H_S3_P2 3 4H

48H_S1_P1 1 48H

48H_S1_P2 1 48H

48H_S2_P1 2 48H

48H_S2_P2 2 48H

48H_S3_P1 3 48H

48H_S3_P2 3 48H

AT_S1_P1 1 ACT

AT_S1_P2 1 ACT

AT_S2_P1 2 ACT

AT_S2_P2 2 ACT

AT_S3_P1 3 ACT

AT_S3_P2 3 ACT

TP_2011_S1_P1 1 TOP

TP_2011_S1_P2 1 TOP

TP_2011_S2_P1 2 TOP

TP_2011_S2_P2 2 TOP

TP_2011_S3_P1 3 TOP

TP_2011_S3_P2 3 TOP

TP_2015_S1_P1 4 TOP

TP_2015_S1_P2 4 TOP

TP_2015_S2_P1 5 TOP

TP_2015_S2_P2 5 TOP

TP_2015_S3_P1 6 TOP

TP_2015_S3_P2 6 TOP