반응형

DESeq2 에서 multiple condition 수행하기




DESeq2는 여러 컨디션이 있어도 DEG분석이란 결국 pairwise하게 비교해야 하기 때문에 condtion을 구분 해 놓은뒤 결과를 원하는 두 컨디션을 놓고 비교해서 각각 분석해야 한다.

SRA010153 데이터로 전체적인 DEG analysis workflow는 추후에 정리하기로 하고 여기서는 multiple condition의 비교에만 중점을 두고 설명하겠다.

SRA010153은 여러 컨디션이 있지만 Brain과 UHR phi X 데이터만 받은 상태에서 임의로 test 컨디션을 아래처럼 추가했다. (원래는 7개의 Brain과 7개의 UHR이다.)

sampleCondition<-c('Brain','Brain','Brain','Brain','Brain','Brain','test', 'UHR','UHR','UHR','UHR','test','UHR','test')

sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles, condition=sampleCondition)    

                                                                   

ddsHTSeq<-DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory=directory, design=~condition)

write.csv(assay(ddsHTSeq),quote=F,file="raw_count.csv")                                                                               

                                                

colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c('Brain','UHR','test'))

                                                                                        

cts <- counts(ddsHTSeq)                                                                 

geoMeans <- apply(cts, 1, function(row) if (all(row == 0)) 0 else exp(mean(log(row[row != 0]))))

ddsHTSeq <- estimateSizeFactors(ddsHTSeq, geoMeans=geoMeans)                            

dds<-DESeq(ddsHTSeq)


res <- results(dds)


쭉 진행하면 나오는 결과는 test와 brain간의 비교이다. DESeq2는 default로 가장 처음과 마지막에 있는 컨디션을 비교하기 때문이다.


이를 아래처럼 수정하면 test와 UHR간의 비교 결과를 보여준다.


results(dds,contrast=c('condition','test','UHR'))


2, 3번째 값을 수정하면 그에 맞게 바뀌는데 당연히 이미 설정해놓은 컨디션 값을 지정해주어야하며


순서를 바꾸면 log2FoldChange값이 +에서 -로 바뀌게 된다.


> results(dds,contrast=c('condition','UHR','test'))

log2 fold change (MLE): condition UHR vs test 

Wald test p-value: condition UHR vs test 

DataFrame with 60461 rows and 6 columns

                          baseMean     log2FoldChange             lfcSE

                         <numeric>          <numeric>         <numeric>

ENSG00000000003   92.1871316897525  0.479423175690374  0.43657650183317

ENSG00000000005   2.33464385903914  0.243364522893621 0.807778836598724


> results(dds,contrast=c('condition','test','UHR'))

log2 fold change (MLE): condition test vs UHR 

Wald test p-value: condition test vs UHR 

DataFrame with 60461 rows and 6 columns

                          baseMean      log2FoldChange             lfcSE

                         <numeric>           <numeric>         <numeric>

ENSG00000000003   92.1871316897525  -0.479423175690374  0.43657650183317

ENSG00000000005   2.33464385903914  -0.243364522893621 0.807778836598724


log2FC값 외에 변경되는 값은 없다.





반응형

+ Recent posts