I am interested in calculating dominance coefficients of gene expression. For this we have RNAseq data consisting of two homozygous lines (i.e. Genotype
BB, respectively) and of the heterozygous crosses between them (i.e. Genotype
BA, respectively). We have data for both sexes in each group and I am using DESeq2 to analyse the data using the following design
design = ~Sex+Genotype. For now I will ignore the
Sex component as it is not relevant for my question.
In a first step I am identifying genes that show different expression in
BB and I can easily do so using the suggested DE steps (thanks to the great DESeq2 manual/tutorial :)).
In a second step I am interested whether there are any heterozygotes (i.e.
AB) that show non-additive gene expression, in other words are there cases where
AB gene expression is different from the mean gene expression between
BB. I can calculated this via the dominance coefficient as 1-(
BB) which would be 1 if
A gene expression is dominant over
B, 0.5 if they are additive and 0 if
A is recessive to
B. This overall seems to work well and I get a nice distribution around 0.5 meaning that the majority of genes show additive gene expression.
Is there a way to use DESeq2 to test if a gene shows significant dominance or recessivity? I.e. a way to do DE between
AB and the mid-gene expression between
My solution so far is to manually create all homozygote mid-gene expressions (that I call
AABB) where I take the mean between any
BB pair. I.e. I have 3 ´AA´ samples and 2
BB samples, I would create 6
AA1BB1, AA1BB2, AA2,BB1, AA2BB2, AA3BB1, AA3BB2. I do the same thing for the heterozygotes where I similarly pair up
BA to get
ABBA and then do the DE contrast
ABBA. This sort of works as it seems to gives me reasonable results (as far as I can judge) but I am not happy with my solution because I am concerned that this is not the best way or potentially a really bad way of doing this as it may give DESeq2 false power due to potential pseudo replication?
Any help would be greatly appreciated