DESeq雙因素差異分析

單因素差異分析

對于DESeq來說熏挎，常用的單因素的差異分析已經(jīng)比較多了赃蛛，



一般來說，對于基因的表達矩陣富俄，放在前面的為對照組缕题，后面的為處理組截歉，所以結(jié)果為后面的比較前面的
并且對于control組或者treat組這樣的組來說，組內(nèi)各樣本的表達譜姑且認為是相同的

雙因素差異分析

那么烟零，對于雙因素差異分析



所謂雙因素瘪松，是指一個為主因素，另一個為次要因素锨阿，我們除了考慮treat與control的差異以外凉逛，分別在這兩個組內(nèi)，不同genotype也會造成彼此之間的差異

雙因素差異分析的交互項：
如果不考慮交互項群井，那么意味著在 genotype I 状飞，genotype II 和 genotype III 中 treat 比較于 control 的差異值都是相等的；而事實并非如此，在 genotype I 诬辈，genotype II 和 genotype III 中 treat 比較于 control 的差異值很可能是不一樣的酵使，因此需要引入交互項來評價在 genotype I ，genotype II 和 genotype III 中 treat 比較于 control 差異值的不同


交互項的意義是當控制第二個變量時焙糟，我們討論 control 與treat 的差異口渔，那么效應值為 a2，這意味著在第二個變量 genotype I 穿撮，genotype II 和 genotype III 中 control 與treat 的差異都是 a2缺脉。可是不是在第二個變量 genotype I 悦穿，genotype II 和 genotype III 類中的差異都是 a2攻礼？需要引入交互項來區(qū)分這種差異
因此引入交互項來區(qū)分在第二個變量 genotype I ，genotype II 和 genotype III 中 control 與treat 的差異栗柒，如上圖 y2 和 y3礁扮，同樣描述 control 與treat 的差異（效應值為 a2）但在第二個變量 genotype II 和 genotype III 中是不同的，在 genotype II 類中 control 與 treat 的差異為 a2+d5瞬沦，在 genotype III 中 control 與 treat 的差異為a2+d6


交互項顯著代表 control 與 treat 的差異在不同的 genotype 是不同的太伊；否則相同

雙因素差異分析結(jié)構(gòu)


這就是各個因素之間的結(jié)構(gòu)和內(nèi)在聯(lián)系，一般我們把設計矩陣的前面那幾列作為對照

雙因素實例

若不帶交互項
dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$genotype <- factor(rep(rep(c("I","II","III"),each=3),2))
#我們把主要因素condition寫在后面逛钻，genotype為次要因素寫在前面僚焦，“：”表示交互項
design(dds) <- ~ genotype + condition
dds <- DESeq(dds)
resultsNames(dds)

我們可以發(fā)現(xiàn)，若沒有交互項曙痘，各組的比較僅停留在單因素兩兩進行比較芳悲。
設計矩陣的因子水平如下：


設計矩陣的因子水平

第一種情況：

dat = counts(dds, normalized=TRUE)
dds_result = data.frame(results(dds, contrast=c("condition","B","A")))

這樣的差異表達方式只考慮condition B和A的差異，而不區(qū)分是哪一種genotype
我們挑選前三個基因手動計算它們的log2(FC)值：

# dat 的 10-18列包括了所有 genotype I-III 的 conditionB 的 sample
# dat 的 1-9列包括了所有 genotype I-III 的 conditionA 的 sample

## Gene 1
log2(mean(dat[1,10:18])/mean(dat[1,1:9]))
[1] -0.5054413

## Gene 2
log2(mean(dat[2,10:18])/mean(dat[2,1:9]))
[1] -0.09992514

## Gene 3
log2(mean(dat[3,10:18])/mean(dat[3,1:9]))
[1] -0.2876421

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.495586066
  2. Gene 2: log2(FC) 為 -0.103780171
  3. Gene 3: log2(FC) 為 -0.322288959

第二種情況：

dat = counts(dds, normalized=TRUE)
dds_result = data.frame(results(dds, contrast=c("genotype","II","I")))

這樣的差異表達方式只考慮genotype II 和 I 之間的差異屡江，而不區(qū)分是哪一種condition
我們挑選前三個基因手動計算它們的log2(FC)值：

# dat 的 4-6列和13-15列包括了所有 condition 的 genotype II 的 sample
# dat 的 1-3列和10-12列包括了所有 condition 的 genotype I 的 sample

## Gene 1
log2(mean(dat[1,c(4:6,13:15)])/mean(dat[1,c(1:3,10:12)]))
[1] -0.3342552

# Gene 2
log2(mean(dat[2,c(4:6,13:15)])/mean(dat[2,c(1:3,10:12)]))
[1] 0.4529269

# Gene 3
log2(mean(dat[3,c(4:6,13:15)])/mean(dat[3,c(1:3,10:12)]))
[1] -0.1632169

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.3179995795
  2. Gene 2: log2(FC) 為 0.4551878042
  3. Gene 3: log2(FC) 為 -0.2035258259

第三種情況：

dat = counts(dds, normalized=TRUE)
dds_result = data.frame(results(dds, contrast=c("genotype","III","I")))

這樣的差異表達方式只考慮genotype III 和 I 之間的差異芭概，而不區(qū)分是哪一種condition
我們挑選前三個基因手動計算它們的log2(FC)值：

# dat 的 7-9列和16-18列包括了所有 condition 的 genotype III 的 sample
# dat 的 1-3列和10-12列包括了所有 condition 的 genotype I 的 sample

## Gene 1
log2(mean(dat[1,c(7:9,16:18)])/mean(dat[1,c(1:3,10:12)]))
[1] -0.5673898

# Gene 2
log2(mean(dat[2,c(7:9,16:18)])/mean(dat[2,c(1:3,10:12)]))
[1] 0.7163037

# Gene 3
log2(mean(dat[3,c(7:9,16:18)])/mean(dat[3,c(1:3,10:12)]))
[1] -0.2464887

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.5515813451
  2. Gene 2: log2(FC) 為 0.7185941762
  3. Gene 3: log2(FC) 為 -0.3108995692

第四種情況：

dat = counts(dds, normalized=TRUE)
dds_result = data.frame(results(dds, contrast=c("genotype","III","II")))

這樣的差異表達方式只考慮genotype III 和 II 之間的差異，而不區(qū)分是哪一種condition
我們挑選前三個基因手動計算它們的log2(FC)值：

# dat 的 7-9列和16-18列包括了所有 condition 的 genotype III 的 sample
# dat 的 4-6列和13-15列包括了所有 condition 的 genotype II 的 sample

## Gene 1
log2(mean(dat[1,c(7:9,16:18)])/mean(dat[1,c(4:6,13:15)]))
[1] -0.2331346

# Gene 2
log2(mean(dat[2,c(7:9,16:18)])/mean(dat[2,c(4:6,13:15)]))
[1] 0.2633768

# Gene 3
log2(mean(dat[3,c(7:9,16:18)])/mean(dat[3,c(4:6,13:15)]))
[1] -0.08327181

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.23358177
  2. Gene 2: log2(FC) 為 0.26340637
  3. Gene 3: log2(FC) 為 -0.10737374
若帶有交互項

我們用幫助文檔中的例子作為說明

dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$genotype <- factor(rep(rep(c("I","II","III"),each=3),2))
#我們把主要因素condition寫在后面惩嘉，genotype為次要因素寫在前面罢洲，“：”表示交互項
design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)

如果帶有交互項，那么我們就可以比較雙因素共同作用下兩組的差異了文黎；也可以限定一個因素惹苗，比較某單因素作用下兩組的差異

注意事項

DESeq2的建模原理是利用廣義線性模型去擬合高維數(shù)據(jù)，然后再做差異比較
值得注意的是耸峭，不帶有交互項（design(dds) <- ~ genotype + condition）對高維數(shù)據(jù)的擬合的方程是線性的（兩個維度桩蓉，其中一個維度是genotype，另一個維度是condition）劳闹；帶有交互項（design(dds) <- ~ genotype + condition + genotype:condition）對高維數(shù)據(jù)的擬合的方程是非線性的（兩個維度院究，其中一個維度是genotype洽瞬，另一個維度是condition ，再帶有一個非線性交互項 genotype:condition）业汰；因此兩個模型計算出來的基因表達估計值是不一樣的伙窃，因此計算出來的log2(FC)是不一樣的，但是不會影響上下調(diào)的方向

上面是沒有交互項的样漆；下面是有交互項的

那么DESeq的原理是哪一個sample放前面为障，哪一個就為對照組，在該例子里面放祟，對于condition鳍怨，condition_A為對照，對于genotype跪妥，genotype I為對照組
當然我們也可以自行設定：

# to compare A vs B, make B the reference level,
# and select the last coefficient
condition <- relevel(condition, "B")
  1. genotype_II_vs_I 表示的是在condition_A的條件下鞋喇，genotype II 與genotype I 相比較下的差異
  2. genotype_III_vs_I 表示的是在condition_A的條件下，genotype III 與genotype I 相比較下的差異
  3. condition_B_vs_A 表示對于genotype I 來說骗奖，在condition_B與condition_A相比較下的差異
  4. genotypeII.conditionB 表示雙因素交互項确徙，目的是檢測在不同genotype中醒串，condition的效果是否會不同执桌，即先計算condition_B的genotype II 與condition_A的genotype II 相比較下的差異和condition_B的genotype I 與condition_A的genotype I 的差異，兩者再進行差異比較（前者比后者）
  5. genotypeIII.conditionB 表示雙因素交互項芜赌，目的是檢測在不同genotype中仰挣，condition的效果是否會不同，即先計算condition_B的genotype III 與condition_A的genotype III 相比較下的差異和condition_B的genotype I 與condition_A的genotype I 的差異缠沈，兩者再進行差異比較（前者比后者）

接下來我們介紹下參數(shù) contrast：


這個參數(shù)目的是指定要從對象中提取哪些比較以構(gòu)建結(jié)果表
簡而言之膘壶，這個參數(shù)控制著提取誰與誰比較的結(jié)果

第一種情況：
那么：

dat = counts(dds, normalized=TRUE)
res = results(dds, list( c("genotype_II_vs_I")))

這個表示對于condition_A 來說，在genotype II 與genotype I 相比較下的差異
FC＝(genotype II in condition_A) / (genotype I in condition_A)
我們用代碼計算下：

# dat 的 4-6列代表 genotype II in condition_A 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2(mean(dat[1,4:6])/mean(dat[1,1:3]))
[1] 0.3768983

## Gene 2
log2(mean(dat[2,4:6])/mean(dat[2,1:3]))
[1] 0.3134871

## Gene 3
log2(mean(dat[3,4:6])/mean(dat[3,1:3]))
[1] -0.5114681

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 0.367333
  2. Gene 2: log2(FC) 為 0.313593
  3. Gene 3: log2(FC) 為 -0.506202

第二種情況：
那么：

dat = counts(dds, normalized=TRUE)
res = results(dds, list( c("genotype_III_vs_I")))

這個表示對于condition_A 來說洲愤，在genotype III 與genotype I 相比較下的差異
FC＝(genotype III in condition_A) / (genotype I in condition_A)
我們用代碼計算下：

# dat 的 7-9列代表 genotype III in condition_A 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2(mean(dat[1,7:9])/mean(dat[1,1:3]))
[1] -0.3639434

## Gene 2
log2(mean(dat[2,7:9])/mean(dat[2,1:3]))
[1] -0.03030897

## Gene 3
log2(mean(dat[3,7:9])/mean(dat[3,1:3]))
[1] -0.6220792

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.364814
  2. Gene 2: log2(FC) 為 -0.029441
  3. Gene 3: log2(FC) 為 -0.610125

第三種情況：
那么：

dat = counts(dds, normalized=TRUE)
# the condition effect for genotype I (the main effect)
results(dds, contrast=c("condition","B","A"))

這個表示對于genotype I 來說颓芭，在condition_B與condition_A相比較下的差異
FC＝(genotype I in condition_B) / (genotype I in condition_A)
我們用代碼計算下：

# dat 的 10-12列代表 genotype I in condition_B 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2(mean(dat[1,10:12])/mean(dat[1,1:3]))
[1] -0.6302513

## Gene 2
log2(mean(dat[2,10:12])/mean(dat[2,1:3]))
[1] -0.1265862

## Gene 3
log2(mean(dat[3,10:12])/mean(dat[3,1:3]))
[1] 0.3226792

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.633304150
  2. Gene 2: log2(FC) 為 -0.129958993
  3. Gene 3: log2(FC) 為 0.321771376

第四種情況：

dat = counts(dds, normalized=TRUE)
# the condition effect for genotype III.
# this is the main effect *plus* the interaction term
# (the extra condition effect in genotype III compared to genotype I).
results(dds, list( c("condition_B_vs_A","genotypeIII.conditionB")))

表示主要考察在genotype III 這個因素下，condition_B 與condition_A 相比帶來的差異柬赐。
也就是說亡问，condition_B的genotype III 與 condition_A的genotype III 相比較所帶來的差異
這里的extra condition effect應該理解為在genotype III 的條件下，某兩個condition相比較下的差異
FC ＝ (genotype III in condition_B) / (genotype III in condition_A)

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 7-9列代表 genotype III in condition_A 的 sample

## Gene 1
log2(mean(dat[1,16:18])/mean(dat[1,7:9]))
[1] -0.4690318

## Gene 2
log2(mean(dat[2,16:18])/mean(dat[2,7:9]))
[1] -0.0880866

## Gene 3
log2(mean(dat[3,16:18])/mean(dat[3,7:9]))
[1] -0.8147309

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.46013841
  2. Gene 2: log2(FC) 為 -0.08729557
  3. Gene 3: log2(FC) 為 -0.81238535

第五種情況：

dat = counts(dds, normalized=TRUE)
# the condition effect for genotype II.
# this is the main effect *plus* the interaction term
# (the extra condition effect in genotype II compared to genotype I).
results(dds, list( c("condition_B_vs_A","genotypeII.conditionB")))

表示主要考察在genotype II 這個因素下肛宋，condition_B 與condition_A 相比帶來的差異州藕。
也就是說，condition_B的genotype II 與 condition_A的genotype II 相比較所帶來的差異
這里的extra condition effect應該理解為在genotype II 的條件下酝陈，某兩個condition相比較下的差異
FC ＝ (genotype II in condition_B) / (genotype II in condition_A)

# dat 的 13-15列代表 genotype II in condition_B 的 sample
# dat 的 4-6列代表 genotype II in condition_A 的 sample

## Gene 1
log2(mean(dat[1,13:15])/mean(dat[1,4:6]))
[1] 1.375654

## Gene 2
log2(mean(dat[2,13:15])/mean(dat[2,4:6]))
[1] -0.3227521

## Gene 3
log2(mean(dat[3,13:15])/mean(dat[3,4:6]))
[1] 0.08036955

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 1.3656608
  2. Gene 2: log2(FC) 為 -0.3305137
  3. Gene 3: log2(FC) 為 0.0780536

第六種情況：

dat = counts(dds, normalized=TRUE)
# the interaction term for condition effect in genotype III vs genotype I.
# this tests if the condition effect is different in III compared to I
results(dds, name="genotypeIII.conditionB")

這一項是交互項床玻，即目的是檢測在不同genotype中，condition的效果是否會不同沉帮，即先計算condition_B的genotype III 與condition_A的genotype III 相比較下的差異和condition_B的genotype I 與condition_A的genotype I 的差異锈死，兩者再進行差異比較（前者比后者）贫堰，英文的解釋為：This tests if the condition effect is different in III compared to I:

FC ＝ (genotype III in condition_B / genotype III in condition_A) / (genotype I in condition_B / genotype I in condition_A)

這個 FC 主要考察的是 condition 的改變對 genotype III 產(chǎn)生的變化倍數(shù)是否和 condition 的改變對 genotype I 產(chǎn)生的變化倍數(shù)相同；如果相同 FC的值為1待牵；如果前者變化倍數(shù)比后者小严嗜，那么 FC 介于0到1，如果前者變化倍數(shù)比后者大洲敢，那么 FC 大于1
代碼如下：

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 7-9列代表 genotype III in condition_A 的 sample
# dat 的 10-12列代表 genotype I in condition_B 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2((mean(dat[1,16:18])/mean(dat[1,7:9]))/(mean(dat[1,10:12])/mean(dat[1,1:3])))
[1] 0.1612196

## Gene 2
log2((mean(dat[2,16:18])/mean(dat[2,7:9]))/(mean(dat[2,10:12])/mean(dat[2,1:3])))
[1] 0.03849963

## Gene 3
log2((mean(dat[3,16:18])/mean(dat[3,7:9]))/(mean(dat[3,10:12])/mean(dat[3,1:3])))
[1] -1.13741

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 0.17316574
  2. Gene 2: log2(FC) 為 0.04266342
  3. Gene 3: log2(FC) 為 -1.13415673

這里在插播一下name這個參數(shù)：


通常只提取一項的話漫玄，用 name 比較好

第七種情況：

dat = counts(dds, normalized=TRUE)
# the interaction term for condition effect in genotype III vs genotype II.
# this tests if the condition effect is different in III compared to II
results(dds, contrast=list("genotypeIII.conditionB", "genotypeII.conditionB"))

這個是兩個交互項的比較，其英文解釋為：This tests if the condition effect is different in III compared to II
目的是檢測在不同genotype中压彭，condition的效果是否會不同睦优，即先計算condition_B的genotype III 與condition_A的genotype III 相比較下的差異和condition_B的genotype II 與condition_A的genotype II 的差異，兩者再進行差異比較（前者比后者）
FC ＝ (genotype III in condition_B / genotype III in condition_A) / (genotype II in condition_B / genotype II in condition_A)

這個 FC 主要考察的是 condition 的改變對 genotype III 產(chǎn)生的變化倍數(shù)是否和 condition 的改變對 genotype II 產(chǎn)生的變化倍數(shù)相同壮不；如果相同 FC的值為1汗盘；如果前者變化倍數(shù)比后者小，那么 FC 介于0到1询一，如果前者變化倍數(shù)比后者大隐孽，那么 FC 大于1

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 7-9列代表 genotype III in condition_A 的 sample
# dat 的 13-15列代表 genotype II in condition_B 的 sample
# dat 的 4-6列代表 genotype II in condition_A 的 sample

## Gene 1
log2((mean(dat[1,16:18])/mean(dat[1,7:9]))/(mean(dat[1,13:15])/mean(dat[1,4:6])))
[1] -0.08777481

## Gene 2
log2((mean(dat[2,16:18])/mean(dat[2,7:9]))/(mean(dat[2,13:15])/mean(dat[2,4:6])))
[1] 0.006579703

## Gene 3
log2((mean(dat[3,16:18])/mean(dat[3,7:9]))/(mean(dat[3,13:15])/mean(dat[3,4:6])))
[1] -0.3210265

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.072876535
  2. Gene 2: log2(FC) 為 0.007112289
  3. Gene 3: log2(FC) 為 -0.321480534

第八種情況：

dat = counts(dds, normalized=TRUE)
res = results(dds, list( c("genotype_III_vs_I","genotypeIII.condition_B")))

這個比較的是在condition_B的條件下，genotype III 和genotype I 相比較下的差異
即condition_B的genotype III 與condition_B的genotype I 相比較下的差異
FC = (genotype III in condition_B) / (genotype I in condition_B)
手動計算代碼：

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 10-12列代表 genotype I in condition_B 的 sample

## Gene 1
log2(mean(dat[1,16:18])/mean(dat[1,10:12]))
[1] -0.4716114

## Gene 2
log2(mean(dat[2,16:18])/mean(dat[2,10:12]))
[1] 0.7362693

## Gene 3
log2(mean(dat[3,16:18])/mean(dat[3,10:12]))
[1] -0.8629554

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.46202782
  2. Gene 2: log2(FC) 為 0.73984754
  3. Gene 3: log2(FC) 為 -0.86640184

第九種情況：

dat = counts(dds, normalized=TRUE)
res = results(dds, list( c("genotype_II_vs_I","genotypeII.conditionB")))

這個比較的是在condition_B的條件下健蕊，genotype II 和genotype I 相比較下的差異
即condition_B的genotype II 與condition_B的genotype I 相比較下的差異
FC = (genotype II in condition_B) / (genotype I in condition_B)
手動計算代碼：

# dat 的 13-15列代表 genotype II in condition_B 的 sample
# dat 的 10-12列代表 genotype I in condition_B 的 sample

## Gene 1
log2(mean(dat[1,13:15])/mean(dat[1,10:12]))
[1] -0.9762591

## Gene 2
log2(mean(dat[2,13:15])/mean(dat[2,10:12]))
[1] -0.4114295

## Gene 3
log2(mean(dat[3,13:15])/mean(dat[3,10:12]))
[1] 1.126437

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.972029
  2. Gene 2: log2(FC) 為 -0.412843
  3. Gene 3: log2(FC) 為 1.132502

第十種情況：

dat = counts(dds, normalized=TRUE)
# 注意要分為兩個 c() 來寫
res = results(dds, list( c("genotype_III_vs_I","genotypeIII.conditionB"),c("genotype_II_vs_I","genotypeII.conditionB")))

這個比較的是在condition_B的條件下菱阵，genotype III 和genotype II 相比較下的差異
即condition_B的genotype III 與condition_B的genotype II 相比較下的差異
FC = (genotype III in condition_B) / (genotype II in condition_B)
手動計算代碼：

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 13-15列代表 genotype II in condition_B 的 sample

## Gene 1
log2(mean(dat[1,16:18])/mean(dat[1,13:15]))
[1] 0.2984419

## Gene 2
log2(mean(dat[2,16:18])/mean(dat[2,13:15]))
[1] -0.1667908

## Gene 3
log2(mean(dat[3,16:18])/mean(dat[3,13:15]))
[1] -0.02535796

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 0.3033151
  2. Gene 2: log2(FC) 為 -0.1632627
  3. Gene 3: log2(FC) 為 -0.0172830

如果我們想要計算跨condition和跨genotype的差異

第十一種情況：

dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$genotype <- factor(rep(rep(c("I","II","III"),each=3),2))
dds$genotype <- relevel(dds$genotype, "III")
#我們把主要因素condition寫在后面，genotype為次要因素寫在前面缩功，“：”表示交互項
design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)


dat = counts(dds, normalized=TRUE)
res = results(dds, list(c("condition_B_vs_A"),c("genotype_I_vs_III")))

這里計算的是condition_B的genotype III 與condition_A的genotype I 相比較下的差異
FC = (genotype III in condition_B) / (genotype I in condition_A)

手動計算代碼：

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2(mean(dat[1,16:18])/mean(dat[1,1:3]))
[1] -0.2684356

## Gene 2
log2(mean(dat[2,16:18])/mean(dat[2,1:3]))
[1] 0.6341815

## Gene 3
log2(mean(dat[3,16:18])/mean(dat[3,1:3]))
[1] 0.111087

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 -0.267888
  2. Gene 2: log2(FC) 為 0.638842
  3. Gene 3: log2(FC) 為 0.108176

第十二種情況：

dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$genotype <- factor(rep(rep(c("I","II","III"),each=3),2))
dds$genotype <- relevel(dds$genotype, "III")
#我們把主要因素condition寫在后面晴及，genotype為次要因素寫在前面，“：”表示交互項
design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)


dat = counts(dds, normalized=TRUE)
res = results(dds, list(c("condition_B_vs_A"),c("genotype_II_vs_III")))

這里計算的是condition_B的genotype III 與condition_A的genotype II 相比較下的差異
FC = (genotype III in condition_B) / (genotype II in condition_A)

手動計算代碼：

# dat 的 16-18列代表 genotype III in condition_B 的 sample
# dat 的 4-6列代表 genotype II in condition_A 的 sample

## Gene 1
log2(mean(dat[1,16:18])/mean(dat[1,4:6]))
[1] 0.07651671

## Gene 2
log2(mean(dat[2,16:18])/mean(dat[2,4:6]))
[1] 0.966533

## Gene 3
log2(mean(dat[3,16:18])/mean(dat[3,4:6]))
[1] -0.1500791

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 0.0778375
  2. Gene 2: log2(FC) 為 0.9659273
  3. Gene 3: log2(FC) 為 -0.1501982

第十三種情況：

dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$genotype <- factor(rep(rep(c("I","II","III"),each=3),2))
dds$genotype <- relevel(dds$genotype, "II")
#我們把主要因素condition寫在后面嫡锌，genotype為次要因素寫在前面虑稼，“：”表示交互項
design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)


dat = counts(dds, normalized=TRUE)
res = results(dds, list(c("condition_B_vs_A"),c("genotype_I_vs_II")))

這里計算的是condition_B的genotype II 與condition_A的genotype I 相比較下的差異
FC = (genotype II in condition_B) / (genotype I in condition_A)

手動計算代碼：

# dat 的 13-15列代表 genotype II in condition_B 的 sample
# dat 的 1-3列代表 genotype I in condition_A 的 sample

## Gene 1
log2(mean(dat[1,13:15])/mean(dat[1,1:3]))
[1] 1.852217

## Gene 2
log2(mean(dat[2,13:15])/mean(dat[2,1:3]))
[1] 0.9619238

## Gene 3
log2(mean(dat[3,13:15])/mean(dat[3,1:3]))
[1] -0.5547154

而軟件計算的log2(FC)值如下（忽略矯正的情況）：

  1. Gene 1: log2(FC) 為 1.812546
  2. Gene 2: log2(FC) 為 0.959211
  3. Gene 3: log2(FC) 為 -0.553264

關于contract的注意事項

contract這個參數(shù)后面接contract = c("..." , "...")，contract = list(c("..." , "..."))和contract = list(c("..."),c("..."))這幾個所對比的是不一樣的
第一個要求的向量長度為3势木，那么只適用于某一個因素不同的兩組進行對比
例如：


對于上圖的condition_B_vs_A蛛倦，這是一個條件，所以有：

results(dds, contrast=c("condition","B","A"))

那么對于多因素的情況啦桌，比如下面的情況

results(dds, contrast=list( c("condition_B_vs_A","genotypeIII.conditionB") ))

主要考察在genotype III 這個因素下溯壶，condition_B 與condition_A 相比帶來的差異。
也就是說震蒋，condition_B的genotype III 與 condition_A的genotype III 相比較所帶來的差異

那么下面的這種情況：

results(dds, contrast=list( c("genotype_III_vs_I","genotypeIII.conditionB"),c("genotype_II_vs_I","genotypeII.conditionB")))

前一個c()作分子茸塞，后一個c()作分母；
表示的是在condition_B條件下查剖，genotype III 與genotype I 相比較下的差異（分子）和condition_B條件下钾虐，genotype II 與genotype I 相比較下的差異（分母），兩者再進行比較下的差異（再相除）

FC = [(genotype III in condition_B) / (genotype I in condition_B)] / (genotype II in condition_B) / (genotype I in condition_B)] = (genotype III in condition_B) / (genotype II in condition_B)笋庄，相當于約分

參考：《DESeq2說明文檔》
https://www.bioconductor.org/packages/release/bioc/vignettes/DESeq2/inst/doc/DESeq2.html
Bioconductor issue：
https://support.bioconductor.org/p/121674/
主要解釋：
https://rstudio-pubs-static.s3.amazonaws.com/329027_593046fb6d7a427da6b2c538caf601e1.html

最后編輯于 ：
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
  • 人面猴
    序言：七十年代末效扫，一起剝皮案震驚了整個濱河市倔监，隨后出現(xiàn)的幾起案子，更是在濱河造成了極大的恐慌菌仁，老刑警劉巖浩习，帶你破解...
    沈念sama閱讀 206,126評論 6贊 481
  • 死咒
    序言：濱河連續(xù)發(fā)生了三起死亡事件，死亡現(xiàn)場離奇詭異济丘，居然都是意外死亡谱秽，警方通過查閱死者的電腦和手機，發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 88,254評論 2贊 382
  • 救了他兩次的神仙讓他今天三更去死
    文/潘曉璐 我一進店門摹迷，熙熙樓的掌柜王于貴愁眉苦臉地迎上來疟赊，“玉大人，你說我怎么就攤上這事峡碉〗矗” “怎么了？”我有些...
    開封第一講書人閱讀 152,445評論 0贊 341
  • 道士緝兇錄：失蹤的賣姜人
    文/不壞的土叔 我叫張陵鲫寄，是天一觀的道長吉执。 經(jīng)常有香客問我，道長地来，這世上最難降的妖魔是什么戳玫？ 我笑而不...
    開封第一講書人閱讀 55,185評論 1贊 278
  • ?港島之戀（遺憾婚禮）
    正文 為了忘掉前任，我火速辦了婚禮靠抑，結(jié)果婚禮上量九，老公的妹妹穿的比我還像新娘适掰。我一直安慰自己颂碧，他們只是感情好，可當我...
    茶點故事閱讀 64,178評論 5贊 371
  • 惡毒庶女頂嫁案：這布局不是一般人想出來的
    文/花漫 我一把揭開白布类浪。 她就那樣靜靜地躺著载城，像睡著了一般。 火紅的嫁衣襯著肌膚如雪费就。 梳的紋絲不亂的頭發(fā)上诉瓦，一...
    開封第一講書人閱讀 48,970評論 1贊 284
  • 城市分裂傳說
    那天，我揣著相機與錄音力细，去河邊找鬼睬澡。 笑死，一個胖子當著我的面吹牛眠蚂，可吹牛的內(nèi)容都是我干的煞聪。 我是一名探鬼主播，決...
    沈念sama閱讀 38,276評論 3贊 399
  • 雙鴛鴦連環(huán)套：你想象不到人心有多黑
    文/蒼蘭香墨 我猛地睜開眼逝慧，長吁一口氣：“原來是場噩夢啊……” “哼昔脯！你這毒婦竟也來了啄糙？” 一聲冷哼從身側(cè)響起，我...
    開封第一講書人閱讀 36,927評論 0贊 259
  • 萬榮殺人案實錄
    序言：老撾萬榮一對情侶失蹤云稚，失蹤者是張志新（化名）和其女友劉穎隧饼，沒想到半個月后，有當?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體静陈，經(jīng)...
    沈念sama閱讀 43,400評論 1贊 300
  • ?護林員之死
    正文 獨居荒郊野嶺守林人離奇死亡燕雁，尸身上長有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點故事閱讀 35,883評論 2贊 323
  • ?白月光啟示錄
    正文 我和宋清朗相戀三年，在試婚紗的時候發(fā)現(xiàn)自己被綠了鲸拥。 大學時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片贵白。...
    茶點故事閱讀 37,997評論 1贊 333
  • 活死人
    序言：一個原本活蹦亂跳的男人離奇死亡，死狀恐怖崩泡，靈堂內(nèi)的尸體忽然破棺而出禁荒，到底是詐尸還是另有隱情，我是刑警寧澤角撞，帶...
    沈念sama閱讀 33,646評論 4贊 322
  • ?日本核電站爆炸內(nèi)幕
    正文 年R本政府宣布呛伴，位于F島的核電站，受9級特大地震影響谒所，放射性物質(zhì)發(fā)生泄漏热康。R本人自食惡果不足惜，卻給世界環(huán)境...
    茶點故事閱讀 39,213評論 3贊 307
  • 男人毒藥：我在死后第九天來索命
    文/蒙蒙 一劣领、第九天 我趴在偏房一處隱蔽的房頂上張望姐军。 院中可真熱鬧，春花似錦尖淘、人聲如沸奕锌。這莊子的主人今日做“春日...
    開封第一講書人閱讀 30,204評論 0贊 19
  • 一樁弒父案村生，背后竟有這般陰謀
    文/蒼蘭香墨 我抬頭看了看天上的太陽惊暴。三九已至，卻和暖如春趁桃，著一層夾襖步出監(jiān)牢的瞬間辽话，已是汗流浹背。 一陣腳步聲響...
    開封第一講書人閱讀 31,423評論 1贊 260
  • 情欲美人皮
    我被黑心中介騙來泰國打工卫病， 沒想到剛下飛機就差點兒被人妖公主榨干…… 1. 我叫王不留油啤，地道東北人。 一個月前我還...
    沈念sama閱讀 45,423評論 2贊 352
  • 代替公主和親
    正文 我出身青樓蟀苛，卻偏偏與公主長得像益咬，于是被迫代替她去往敵國和親。 傳聞我的和親對象是個殘疾皇子屹逛，可洞房花燭夜當晚...
    茶點故事閱讀 42,722評論 2贊 345

推薦閱讀更多精彩內(nèi)容

  • 選擇題部分 1.(),只有在發(fā)生短路事故時或者在負荷電流較大時,變流器中才會有足夠的二次電流作為繼電保護跳閘之用础废。...
    skystarwuwei閱讀 12,750評論 0贊 7
  • 引言指南目的這些指南的主要目的是為一般臨床醫(yī)生提供信息汛骂，指導他們對食管癌的可用診斷/治療策略做出恰當?shù)倪x擇（針對食...
    醫(yī)博云天閱讀 4,262評論 0贊 2
  • ★681★觸覺語顫增強主要見于 E.肺泡內(nèi)炎癥性浸潤 ★682★女，12歲评腺，白血病帘瞭，近日左眼眶出現(xiàn)2cmX2cm無...
    捷龍閱讀 4,506評論 0贊 2
  • 首先在mac中要安裝glew和glfw，這里我安裝這兩個工具使用的是homebrew包管理蒿讥，這東西超級好用蝶念，安裝命...
    yyyqy閱讀 4,847評論 0贊 0
  • 20190121學習筆記： 什么才是受尊敬的企業(yè)？ 一家受到尊敬的公司芋绸， 理所應當為股東創(chuàng)造物質(zhì)財富媒殉， 與此同時，...
    馬唐閱讀 316評論 0贊 0