Waiting
Login processing...

Trial ends in Request Full Access Tell Your Colleague About Jove
Click here for the English version

Cancer Research

RNA 시퀀싱을 위한 3개의 차동 발현 분석 방법: 림마, EdgeR, DESeq2

Published: September 18, 2021 doi: 10.3791/62528
* These authors contributed equally

Summary

RNA 시퀀싱을 위한 차동 발현 분석 방법의 상세한 프로토콜이 제공되었다: 림마, EdgeR, DESeq2.

Abstract

RNA 시퀀싱 (RNA-seq)은 유전 적 변경과 복잡한 생물학적 과정 사이의 관계를 밝힐 수 있으며 종양의 진단, 예후 및 치료에서 큰 가치를 가지고 있기 때문에 전사학에서 가장 널리 사용되는 기술 중 하나입니다. RNA-seq 데이터의 차동 분석은 비정상적인 전사를 식별하는 데 매우 중요하며 림마, EdgeR 및 DESeq2는 차동 분석을 위한 효율적인 도구입니다. 그러나, RNA-seq 차동 분석은 R 언어를 가진 특정 기술과 의학 교육의 교과 과정에서 부족한 적당한 방법을 선택하는 기능이 필요합니다.

본 명세서에서는, 당사는 각각 림마, DESeq2 및 EdgeR을 통해 담랑고아르시노마(CHOL) 및 정상 조직 간의 분화유전자(DEGs)를 식별하고, 그 결과는 화산 플롯및 벤 다이어그램에 도시된다. limma, DESeq2 및 EdgeR의 세 가지 프로토콜은 유사하지만 분석 프로세스 마다 다른 단계가 있습니다. 예를 들어 선형 모델은 림마의 통계에 사용되는 반면 음수 이비알 분포는 edgeR 및 DESeq2에서 사용됩니다. 또한, 정규화된 RNA-seq 카운트 데이터는 EdgeR 및 림마에 필요하지만 DESeq2에는 필요하지 않습니다.

여기서는 림마, EdgeR 및 DESeq2의 세 가지 차동 분석 방법에 대한 자세한 프로토콜을 제공합니다. 세 가지 방법의 결과는 부분적으로 겹칩니다. 세 가지 방법 모두 고유한 장점이 있으며 메서드 선택은 데이터에만 따라 다릅니다.

Introduction

RNA-시퀀싱(RNA-seq)은 많은 장점(예를 들어, 높은 데이터 재현성)을 가진 전사학에서 가장 널리 사용되는 기술 중 하나이며, 복잡한 생물학적 과정의 기능 및 역학에 대한 이해를 크게증가시켰습니다 1,2. 다른 생물학적 맥락에서 비정상적인 전사체의 식별은 또한 분화 유전자 (DEGs)로 알려져 있으며, RNA-seq 분석에서 중요한 단계입니다. RNA-seq는 병인과 관련된 분자 메커니즘 및 생물학적 기능에 대한 깊은 이해를 얻을 수 있게 합니다. 따라서, 차동 분석은종양의진단, 예후 및 치료에 귀중한 것으로 간주되어 왔다3,4,5. 현재, RNA-seq 차동 발현 분석, 특히 림마, DESeq2 및 EdgeR1,6,7을위해 더 많은 오픈 소스 R/바이오 컨덕터 패키지가 개발되었다. 그러나, 차등 분석은 R 언어와 특정 기술과 의료 교육의 교육 과정에서 부족한 적절한 방법을 선택하는 능력이 필요합니다.

본 프로토콜에서, 암 게놈 아틀라스(TCGA)로부터 추출된 담랑고카르시노마(CHOL) RNA-seq 카운트 데이터를 기반으로, 가장 공지된 방법 중 3개(림마8,EdgeR9 및 DESeq210)가각각 R프로그램에 의해 CHOL과 정상 조직 간의 DEGs를 식별하였다. limma, EdgeR 및 DESeq2의 세 가지 프로토콜은 유사하지만 분석 프로세스 마다 다른 단계가 있습니다. 예를 들어, 정규화된 RNA-seq 카운트 데이터는 EdgeR 및 limma8,9에필요하며, DESeq2는 자체 라이브러리 불일치를 사용하여정규화(10)가아닌 데이터를 수정한다. 더욱이, edgeR은 RNA-seq 데이터에 특히 적합하며, 림마는 마이크로어레이 및 RNA-seq에 사용된다. 선형 모델은 LIMma에 의해 채택되어DEGs(12)를평가하고, edgeR의 통계는 경험적 베이추정, 정확한 테스트, 일반화선형 모델 및 준가능성 시험9을포함한 음수 이난 분포를 기반으로 한다.

요약하자면, 우리는 각각 림마, DESeq2 및 EdgeR을 사용하여 RNA-seq 차동 발현 분석의 상세한 프로토콜을 제공합니다. 이 문서를 참조함으로써 사용자는 RNA-seq 차동 분석을 쉽게 수행하고 데이터에 적합한 차동 분석 방법을 선택할 수 있습니다.

Subscription Required. Please recommend JoVE to your librarian.

Protocol

참고: R-스튜디오 프로그램을 열고 R 파일 "DEGs.R"을 로드하면 추가 파일/스크립트에서 파일을 얻을 수 있습니다.

1. 데이터 다운로드 및 사전 처리

  1. 암 게놈 아틀라스(TCGA)에서 콜린고이오사르키노마(CHOL)의 고처리량 시퀀싱(HTSeq) 카운트 데이터를 다운로드한다. 이 단계는 다음 R 코드로 쉽게 수행할 수 있습니다.
    1. R 패키지를 설치하려면 실행을 클릭합니다.
    2. R 패키지를 로드하려면 실행을 클릭합니다.
      if(!requireNamespace("BiocManager", 조용히=TRUE))
      + 설치.패키지 ("바이오 매니저")
      바이오매니저::설치(c("TCGAbiolinks", "요약실험"))
    3. 작업 디렉터리 설정합니다.
      라이브러리 (TCGAbiolinks)
      라이브러리(요약실험)
      setwd("C:/사용자/류시이/데스크탑")
    4. 암 유형을 선택합니다.
      암 <- "TCGA-CHOL"
    5. "GDCquery.R" 파일에서 R 코드를 실행하여 데이터를 다운로드합니다. 파일 "GDCquery.R"은 보충 파일 / 스크립트에서 수집 할 수 있습니다 :
      소스("추가 파일/스크립트/GDCquery.R")
      헤드(cnt)
      ##TCGA-3X-AAVA-01A-11R-A41I-07
      ##ENSG00000000003 4262
      ##ENSG00000000005 1
      ##ENSG00000000419 1254
      ##ENSG00000000457 699
      ##ENSG00000000460 239
      ##ENSG00000000938 334
      참고: 실행 후 CHOLHTSeq 카운트 데이터를 다운로드하여 "cnt"라고 명명하여 행은 앙상블 유전자 ID를 나타내고 컬럼은 샘플 ID를 나타냅니다. 샘플 아이디에서 위치 14-15의 숫자를 주의하십시오. 01에서 09까지의 숫자는 종양을 나타내고 10에서 19까지오는 것은 정상 조직을 나타냅니다.
  2. 앙상블 유전자 ID를 유전자 기호로 변환합니다.
    1. 저장 경로에 따라 어반으로 가져와서 R으로 가져옵니다. 추가 파일(gencode.v22.annotation.gtf)은 보충 파일에서 획득할 수 있습니다.
      gtf_v22 <-rtracklayer::import('보충 파일/gencode.v22.annotation.gtf')
    2. "gtf_v22"에서 R 코드를 실행합니다. 보충 파일 / 스크립트에서 획득 할 수있는 R"파일 :
      출처("보충 파일/스크립트/gtf_v22. R")
    3. 앙상블 유전자 ID를 유전자 기호로 변환하는 기능 "ann"을 적용합니다.
      cnt=ann (cnt,gtf_v22)
  3. 저발현 유전자 필터링
    1. 실행을 클릭하여 R 패키지 "edgeR"을 설치합니다.
      바이오매니저::설치("edgeR")
    2. 달리기를 클릭하여 R 패키지 "edgeR"을 로드합니다.
      라이브러리(에지R)
    3. 다음 R 코드를 실행하여 백만 개당 개(CPM) 값을 가진 유전자를 두 개 이상의 샘플중 하나보다 큰 것으로 유지합니다.
      <-rowSums(cpm(cnt)>1)>=2유지
      cnt <- as.matrix (cnt[keep,])
      참고: 100만(CPM) 값은 읽기 수 대신 사용되며, 다양한 시퀀싱 깊이로 인한 편차를 제거합니다.

2. "림마"를 통한 차분식 분석

  1. 달리기를 클릭하여 R 패키지 "림마"를 설치합니다.
    바이오매니저:설치("림마")
  2. 달리기를 클릭하여 R 패키지 "림마", "edgeR"을 로드합니다.
    도서관(림마)
    라이브러리(에지R)
  3. 다음 R 코드를 실행하여 설계 행렬을 만듭니다.
    그룹 <- substring(colnames(cnt),14,15) # Extract group information
    그룹[그룹 %인% "01"] <- "Cancer" # set '01' as tumor tissue
    그룹[그룹%인% "11"] <- "Normal" # set '11' as normal tissue
    그룹 <- factor (group, levels = c("Normal","Cancer"))
    1. 디자인 행렬을 만듭니다.
      설계 <-model.matrix(~그룹)
      <-콜네임(cnt)
    2. DGEList 개체를 만듭니다.
      dge <-DGEList(카운트 = cnt, 그룹 = 그룹)
    3. 데이터를 정규화합니다.
      dge <- calcNormFactors (dge, 방법 = "TMM")
    4. 다음 R 코드를 실행하여 림마 추세 메서드 기반 차동 식 해석을 수행합니다.
      dge
      #클래스의 #An 객체 "DGEList"
      ##$counts
      ##TCGA-3X-AAVA-01A-11R-A41I-07
      ##TSPAN6 4262
      ##DPM1 1254
      ##SCYL3 699
      ##C1orf112 239
      ##FGR 334
    5. CPM 값을 계산합니다.
      로그지 <-cpm(dge, log=TRUE, 이전.count=3)
    6. 선형 모델에 맞게 실행을 클릭하여 데이터를 예측하거나 변수 간의 관계를 추론합니다.
      핏 <-lmFit (로지, 디자인)
    7. 베이지안을 기준으로 T 값, F 값 및 로그 배당률을 계산합니다.
      <- 이베이즈(적합, 트렌드=TRUE)
    8. 결과 테이블을 추출합니다.
      res_limma<-as.data.frame(topTable(적합,n=Inf))

      헤드(res_limma)
      ## 로그FC AveExpr t P.Value adj. P.Val B
      ##RP11-252E2.2 -4.899493 -2.488589 -20.88052 2.386656e-25 4.931786e-21 47.28823
      ##BX842568.1 -4.347930 -2.595205 -20.14532 1.082759e-24 1.118706e-20 45.83656
      ##CTC-537E7.3 -5.154894 -2.143292 -19.59571 3.452354e-24 2.216114e-20 44.72001
      ##RP11-468N14.3 -6.532259 -2.029714 -19.49409 4.289807e-24 2.216114e-20 44.51056
      ##AP006216.5 -4.507051 -2.670915 -19.25649 7.153356e-24 2.956339e-20 44.01704
      ##RP11-669E14.4 -4.107204 -2.828311 -18.93246 1.448209e-23 4.987633e-20 43.33543
      #The 차분발 분석 결과는 유전자 id, log2 접이식 변화값(logFC), 실험에서 유전자의 평균 로그2 발현 수준(AveExpr), 수정된 t 통계(t), p.값(P.Value), 거짓 발견률(FDR) 수정된 p값(adj)을 포함하는 "res_limma"에 저장된다. P.Val) 및 차별화된 발현 유전자의 로그 확률(B)
      참고: "edgeR"의 함수 "calcNormFactors()"는 샘플 준비 또는 라이브러리 구성 및 시퀀싱으로 인한 영향을 제거하기 위해 데이터를 정규화하는 데 사용되었습니다. 설계 매트릭스의 시공에서, 매트릭스의 ID를 샘플링하기 위해 실험 설계(예를 들어, 조직 유형: 정상 또는 종양 조직)를 일치시켜야 한다. 림마 추세는 시퀀싱 깊이가 동일하는 데이터에 적합하며 림마 붐은 적합합니다: (i) 샘플 라이브러리 크기가 다른 경우; (ii) TMM에 의해 정규화되지 않은 데이터; (iii) 데이터에는 많은 "노이즈"가 있습니다. 양성 로그FC는 유전자가 실험에서 업 규제된다는 것을 의미하며, 음수수는 유전자가 다운 규제된다는 것을 의미합니다.
    9. DEGs를 식별합니다.
      res_limma$sig <- as.factor(
      ifelse (res_limma $ adj. P.Val < 0.05 및 복근(res_limma$logFC) > 2,
      ifelse (res_limma $logFC > 2,''up','down'),'not')) # adj.p 값 < 0.05 및 |log2FC| >= 2는 DEG를 식별하는 임계값입니다.
      요약(res_limma$sig)
      ##down
      ##1880 ​17341 1443
    10. 결과 테이블을 파일에 출력합니다.
      쓰기.csv(res_limma, 파일 = 'result_limma.csv')
    11. 달리기를 클릭하여 R 패키지 "ggplot2"를 설치합니다.
      설치.패키지("ggplot2")
    12. 달리기를 클릭하여 R 패키지 "ggplot2"를 로드합니다.
      라이브러리(ggplot2)
    13. "화산"에서 R 코드를 실행합니다. R"화산 플롯을 만들 수 있습니다. 파일 "화산. R"은 보충 파일에서 획득할 수 있습니다.
      출처("보충 파일/스크립트/화산. R")
      화산 (res_limma"로그FC","adj. P.Val",2,0.05)
      참고: 유전자는 log2FC 및 adj-p 값에 따라 서로 다른 위치에 매핑할 수 있으며, 상향 조절 된 DEG는 빨간색으로 칠하고, 하향 조절 된 DEG는 녹색으로 칠해로 표시됩니다.
    14. 화산 플롯을 저장하려면 내보내기를 클릭합니다.
      참고: 화산 플롯은 다양한 형식(예: PDF, TIFF, PNG, JPEG 형식)으로 생성되고 다운로드할 수 있습니다. 유전자는 그들의 log2FC 및 adj p 값에 따라 다른 위치로 매핑될 수 있고, 상향 조절 DEGs (log2FC > 2, adj p < 0.05)는 빨간색으로 칠하게되며, 다운 규제 DEGs (log2FC < -2, adj p < 0.05)는 녹색으로 착색되고 비 DEG는 회색으로 채색됩니다.

3. "edgeR"을 통한 차동 식 분석

  1. 달리기를 클릭하여 R 패키지 "edgeR"을 로드합니다.
    라이브러리(에지R)
  2. 다음 R 코드를 실행하여 디자인 행렬을 만듭니다.
    그룹 < 서브스트링(콜네임(cnt),14,15)
    그룹[그룹%인% "01"] <- "암"
    그룹[그룹%인% "11"] <- "보통"
    그룹 =인자(그룹, 레벨 = c("정상","암"))
    설계 < 모델.매트릭스(~그룹)
    로네임(디자인) = 콜네임(cnt)
  3. 실행을 클릭하여 DGEList 개체를 만듭니다.
    dge <- DGEList (카운트 = cnt)
  4. 데이터를 정규화합니다.
    dge <- calcNormFactors (dge, 방법 = "TMM")
  5. 실행을 클릭하여 유전자 발현 값의 분산을 추정합니다.
    dge <-추정Disp (dge, 디자인, 견고한 = T)
  6. 실행을 클릭하여 모델을 사용하여 데이터를 계산합니다.
    핏 <-glmQLFit (dge, 디자인)
  7. 통계 테스트를 수행합니다.
    핏 <-glmQLFTest(적합)
  8. 결과 테이블을 추출합니다. 결과는 로그 폴드 변경 값, 로그 CPM, F, p 값 및 FDR 수정 p 값을 포함하는 "res_edgeR"에 저장됩니다.
    res_edgeR=as.data.frame(topTags(적합, n=Inf))
    헤드(res_edgeR)
    ## 로그FC 로그CPM F PValue FDR
    ##GCDH -3.299633 5.802700 458.5991 1.441773e-25 2.979280e-21
    ##MSMO1 -3.761400 7.521111 407.0416 1.730539e-24 1.787993e-20R
    ##CL1 -3.829504 5.319641 376.5043 8.652474e-24 5.516791e-20
    ##ADI1 -3.533664 8.211281 372.6671 1.067904e-23 5.516791e-20
    ##KCNN2 -5.583794 3.504017 358.6525 2.342106e-23 9.679455e-20
    ##GLUD1 -3.287447 8.738080 350.0344 3.848408e-23 1.194406e-19
    #The 결과는 로그 폴드 변경 값(logFC), 로그 CPM, F, p 값 및 FDR 수정 p 값을 포함하는 "res_edgeR"에 저장됩니다.
  9. DEGs를 식별합니다.
    res_edgeR$sig = as.factor(
    ifelse(res_edgeR$FDR < 0.05 및 복근(res_edgeR$logFC) > 2,
    ifelse (res_edgeR$logFC > 2,'up','down'),'''하지 않음))
    요약(res_edgeR$sig)
    ##down
    ##1578 15965 3121
  10. 결과 테이블을 파일에 출력합니다.
    쓰기.csv(res_edgeR, 파일 = 'res_edgeR.csv')
  11. 화산 플롯을 만듭니다.
    화산(res_edgeR"로그FC",FDR",2,0.05)
  12. 화산 플롯을 저장하려면 내보내기를 클릭합니다.

4. "DESeq2"를 통한 차분 식 분석

  1. 실행을 클릭하여 R 패키지 "DESeq2"를 설치합니다.
    바이오매니저:설치("DESeq2")
  2. R 패키지 "DESeq2"를 로드하려면 실행을 클릭합니다.
    라이브러리(DESeq2)
  3. 다음 R 코드를 실행하여 그룹화 계수를 결정합니다.
    그룹 < 서브스트링(콜네임(cnt),14,15)
    그룹[그룹%인% "01"] <- "암"
    그룹[그룹%인% "11"] <- "보통"
    그룹 =인자(그룹, 레벨 = c("정상","암"))
  4. DESeqDataSet 개체를 만듭니다.
    dds <-DESeqDataSetFromMatrix (cnt, DataFrame(그룹), 디자인 = ~그룹)
    dds
    ##class: 데세크데이터세트
    ##dim: 20664 45
    ##metadata(1): 버전
    ##assays(1): 카운트
    ##rownames (20664): TSPAN6 DPM1 ... RP11-274B21.13 LINC01144
    ##rowData 이름(0):
    ##colnames (45): TCGA-3X-AAVA-01A-11R-A41I-07 ...
    ##colData 이름(1): 그룹
  5. 분석을 수행합니다.
    dds <- DESeq (dds)
  6. 결과 테이블을 생성합니다.
    res_DESeq2 <-data.frame(결과(dds))

    헤드(res_DESeq2)
    ## 베이스평균 로그2폴드변경 lfcSE 통계 표판 패드
    ##TSPAN6 4704.9243 -0.8204515 0.3371667 -2.433370 1.495899e-02 2.760180e-02
    ##DPM1 1205.9087 -0.3692497 0.1202418 -3.070894 2.134191e-03 4.838281e-03
    ##SCYL3 954.9772 0.2652530 0.2476441 1.071106 2.841218e-01 3.629059e-01
    ##C1orf112 277.7756 0.7536911 0.2518929 2.992109 2.770575e-03 6.101584e-03
    ##FGR 345.8789 -0.6423198 0.3712729 -1.730047 8.362180e-02 1.266833e-01
    ##CFH 27982.3546 -3.8761382 0.5473363 -7.081823 1.422708e-12 1.673241e-11
    참고: 결과는 정규화된 읽기 수(baseMean), 로그 폴드 변경 값(log2FoldChange), 로그 폴드 변경 표준 오류(lfcSE), Wald 통계(통계), 원본 p값(pvalue) 및 수정된 p 값(padj)의 평균을 포함하는 "res_DESeq2"에 저장됩니다.
  7. DEGs를 식별합니다.
    res_DESeq2$sig = as.factor(
    ifelse (res_DESeq2$ padj < 0.05 및 복근 (res_DESeq2 $log2FoldChange) > 2,
    ifelse(res_DESeq2$log2FoldChange > 2,'up','down'),''하지 않음'))
    요약(res_DESeq2$sig)
    ##down
    ##1616 16110 2938
  8. 결과 테이블을 파일에 출력합니다.
    쓰기.csv(res_DESeq2, 파일 = 'res_DESeq2.csv')
  9. 화산 플롯을 만듭니다.
    화산 (res_DESeq2,"log2FoldChange","padj",2,0.05)
  10. 화산 플롯을 저장하려면 내보내기를 클릭합니다.

5. 벤 다이어그램

  1. 실행을 클릭하여 R 패키지 "VennDiagram"를 설치합니다.
    설치.패키지("벤다이어그램")
  2. R 패키지 "VennDiagram"를 로드하려면 실행을 클릭합니다.
    라이브러리(벤다이어그램)
  3. 최대 규제 된 DEGs의 벤 다이어그램을 확인합니다.
    grid.newpage ()
    grid.draw(venn.다이어그램(목록(림마=로네임(res_
    림마[res_limma$sig="up",]),
    edgeR=로네임(res_edgeR[res_edgeR$sig=="up",]),
    DESeq2=로네임(res_DESeq2[res_DESeq2$sig==
    "위로",))
    NULL, 높이 = 3,너비 = 3,단위 = "in",
    col="블랙", lwd=0.3,fill=c("#FF6666","#FFFF00",
    "#993366"),
    알파=c(0.5, 0.5, 0.5), 메인 = "최대 규제 DEGs"))
  4. 내보내기를 클릭하여 Venn 다이어그램을 저장합니다.
  5. 다운 규제 DEGs의 벤 다이어그램을 확인합니다.
    grid.newpage ()
    grid.draw(venn.다이어그램(목록(림마=로네임(res_
    림마[res_limma$sig="down",]),
    edgeR=행이름(res_edgeR[res_edgeR$sig==
    "아래로",]),
    DESeq2=로네임(res_DESeq2[res_DESeq2$sig=="down",])))
    NULL, 높이 = 3,너비 = 3,단위 = "in",
    col="블랙", lwd=0.3,fill=c("#FF6666","#FFFF00",
    "#993366"),
    알파=c(0.5, 0.5, 0.5), 메인 = "다운 규제 DEGs"))
  6. 내보내기를 클릭하여 Venn 다이어그램을 저장합니다.

Subscription Required. Please recommend JoVE to your librarian.

Representative Results

화산 플롯과 Venn 다이어그램이 특히 사용되는 차동 식 분석의 결과를 시각화하는 다양한 접근 방식이 있습니다. 림마는 |로그FC|≥2와 adj를 가진 CHOL과 정상 조직 사이 3323개의 DEG를 확인했습니다. P.Val<0.05는 임계값으로, 그 중 1880은 CHOL 조직에서 하향 조절되었고 1443은 업규제(그림 1a)였다. 한편, EdgeR은 1578개의 다운 규제 DEGs와 3121개의 업 규제 DEGs(그림1b);를 확인했습니다. DESeq2는 1616개의 하향 규제 DEGs와 2938개의 UP-regulated DEGs(그림1c)를 확인했습니다. 이들 세 가지 방법의 결과를 비교하면, 1431개의 상향 조절 된 DEGs 와 1531 개의 다운 규제 DEGs가 겹쳐졌다(그림 2).

Figure 1
그림 1. CHOL과 정상 조직 사이의 차별화된 발현 유전자(DEGs)의 식별. (a-c) 림마, edgeR 및 DESeq2에 의해 획득 된 모든 유전자의 화산 플롯은 각각 접이식 변화 (log2)에 대해 플롯되며, 빨간색 점은 최대 규제 DEGs (조정 된 p 값<0.05 및 로그 | FC|> 2) 및 그린 포인트는 다운 규제 DEGs(조정된 p 값< 0.05 및 로그 | 나타냅니다. FC|< 2). 이 그림의 더 큰 버전을 보려면 여기를 클릭하십시오.

Figure 2
그림 2. Venn 다이어그램은 림마, 에지R 및 DESeq2에서 파생된 결과 중 중복을 보여 준다. 이 그림의 더 큰 버전을 보려면 여기를 클릭하십시오.

보충 파일. 이 파일을 다운로드하려면 여기를 클릭하십시오.

Subscription Required. Please recommend JoVE to your librarian.

Discussion

암에 있는 풍부한 수차성 전사체는 RNA-seq 차동 분석에 의해 쉽게 확인할 수 있습니다5. 그러나, RNA-seq 차동 발현 분석의 적용은 R 언어와 적절한 방법을 선택할 수 있는 특정 능력을 필요로 하기 때문에 종종 제한됩니다. 이 문제를 해결하기 위해, 우리는 RNA-seq 차동 발현 분석을 적용하기위한 세 가지 가장 잘 알려진 방법 (limma, EdgeR 및 DESeq2)에 대한 자세한 소개및 자습서를 제공합니다. 이를 통해 세 가지 방법 모두에 걸쳐 유사점과 차이점을 쉽게 이해하고 개별 데이터에 적합한 방법을 선택할 수 있으며 복잡한 동적 생물학적 과정을 이해할 수 있습니다.

여기서, 우리는 각각 리마, EdgeR 및 DESeq2를 통해 RNA-seq 차동 발현 분석을 위한 상세한 프로토콜을 제시합니다: (i) 데이터의 다운로드 및 사전 처리, (ii-iv) 리마를 통한 차동 발현 분석, (ii-iv) 리마, 에지R 및 DESeq2를 통해 각각, (v) 벤다이어그램을 통해 이들 세 가지 방법의 결과의 비교.

세 가지 방법은 차동 식 분석 의 프로세스 간에 유사하고 다른 단계를 갖습니다. 선형 모델은 마이크로어레이, RNA-seq 및 정량PCR8,13을포함한 모든 유전자 발현 기술에 적용되는 림마의 통계에 사용되며, edgeR 및 DESeq2는 음의 이소성 분포9,10,및 edgeR 및 DESeq2가 RNA-seq 데이터에 적합합니다. 또한, 정규화된 RNA-seq 카운트 데이터는 EdgeR 및 림마에 필요하며, DESeq2는 자체 라이브러리 불일치를 사용하여 정규화 대신 데이터를 수정하고 DESeq2의 데이터는 정수 매트릭스여야 합니다. 정규화 방법에는 TMM(M-값의 트리밍 평균), TMMwsp, RLE(상대 로그 식) 및 상층체를 포함하며, 그 중 TMM은 RNA-seq 데이터에 가장 일반적으로 사용되는 정규화 방법입니다. 세 가지 방법의 결과는 DESeq2 및 EdgeR이 림마보다 더 많은 DEGs를 얻는 것으로 나타났습니다. 이러한 차이의 이유는 edgeR및 DESeq2가 많은 수의 거짓 긍정에 기여하는 음의 이소성 모델을 기반으로 하기 때문입니다. 반대로, 림마붐은 분산 기능만을 사용하며, 리마14,15,16을사용한 선형 모델 분석이 선행된 분산 안정화 변환의 경우와 마찬가지로 과도한 거짓긍정을나타내지 않는다.

세 가지 방법 모두 고유한 장점이 있으며, 선택은 데이터 유형에 만 따라 달라집니다. 예를 들어 마이크로어레이 데이터가 있는 경우 림마를 우선적으로 제공해야 하지만 차세대 시퀀싱 데이터인 경우 DESeq2 및 EdgeR이9,10,17을선호합니다. 요약하자면, R 패키지 림마, 에지R 및 DESeq2를 각각 RNA-seq 차동 발현 분석을 위한 상세한 프로토콜을 제공합니다. 세 가지 메서드의 출력 결과는 부분적으로 겹치며 이러한 차등 메서드는 각각의 장점이 있습니다. 안타깝게도 이 프로토콜은 다른 데이터 유형(예: 마이크로어레이 데이터) 및 메서드(예: EBSeq)18에대한 기술적 세부 정보를 다루지 않습니다.

Subscription Required. Please recommend JoVE to your librarian.

Disclosures

원고는 이전에 출판되지 않았으며 다른 곳에서 출판될 것으로 간주되지 않습니다. 모든 저자는 중요한 지적 콘텐츠에 대한이 원고의 작성에 기여하고 읽고 최종 원고를 승인했다. 우리는 이해 상충이 없다고 선언합니다.

Acknowledgments

이 작품은 중국 국립 자연과학 재단(81860276 보조금)과 국가 핵심 R&D 프로그램의 주요 특별 기금 프로젝트(보조금 2018YFC1003200)의 지원을 받았습니다.

Materials

Name Company Catalog Number Comments
R version 3.6.2 free software
Rstudio free software

DOWNLOAD MATERIALS LIST

References

  1. Tambonis, T., Boareto, M., Leite, V. B. P. Differential Expression Analysis in RNA-seq Data Using a Geometric Approach. Journal of Computational Biology. 25, 1257-1265 (2018).
  2. Wang, Z., Gerstein, M., Snyder, M. RNA-Seq: a revolutionary tool for transcriptomics. Nature Reviews. Genetics. 10, 57-63 (2009).
  3. Anders, S., et al. Count-based differential expression analysis of RNA sequencing data using R and Bioconductor. Nature Protocols. 8, 1765-1786 (2013).
  4. McDermaid, A., Monier, B., Zhao, J., Liu, B., Ma, Q. Interpretation of differential gene expression results of RNA-seq data: review and integration. Briefings in Bioinformatics. 20, 2044-2054 (2019).
  5. Costa-Silva, J., Domingues, D., Lopes, F. M. RNA-Seq differential expression analysis: An extended review and a software tool. PloS One. 12, 0190152 (2017).
  6. Law, C. W., et al. RNA-seq analysis is easy as 1-2-3 with limma, Glimma and edgeR. F1000Research. 5, (2016).
  7. Varet, H., Brillet-Guéguen, L., Coppée, J. Y., Dillies, M. A. SARTools: A DESeq2- and EdgeR-Based R Pipeline for Comprehensive Differential Analysis of RNA-Seq Data. PloS One. 11, 0157022 (2016).
  8. Ritchie, M. E., et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Research. 43, 47 (2015).
  9. Robinson, M. D., McCarthy, D. J., Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 26, Oxford, England. 139-140 (2010).
  10. Love, M. I., Huber, W., Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology. 15, 550 (2014).
  11. Gentleman, R. C., et al. Bioconductor: open software development for computational biology and bioinformatics. Genome Biology. 5, 80 (2004).
  12. Law, C. W., Chen, Y., Shi, W., Smyth, G. K. voom: Precision weights unlock linear model analysis tools for RNA-seq read counts. Genome Biology. 15, 29 (2014).
  13. Smyth, G. K. Linear models and empirical bayes methods for assessing differential expression in microarray experiments. Statistical Applications in Genetics and Molecular Biology. 3, (2004).
  14. Lund, S. P., Nettleton, D., McCarthy, D. J., Smyth, G. K. Detecting differential expression in RNA-sequence data using quasi-likelihood with shrunken dispersion estimates. Statistical Applications in Genetics and Molecular Biology. 11, (2012).
  15. Reeb, P. D., Steibel, J. P. Evaluating statistical analysis models for RNA sequencing experiments. Frontiers in Genetics. 4, 178 (2013).
  16. Rocke, D. M., et al. Excess False Positive Rates in Methods for Differential Gene Expression Analysis using RNA-Seq Data. bioRxiv. , (2015).
  17. Agarwal, A., et al. Comparison and calibration of transcriptome data from RNA-Seq and tiling arrays. BMC genomics. 11, 383 (2010).
  18. Leng, N., et al. EBSeq: an empirical Bayes hierarchical model for inference in RNA-seq experiments. Bioinformatics. 29, Oxford, England. 1035-1043 (2013).

Tags

암 연구 문제 175
RNA 시퀀싱을 위한 3개의 차동 발현 분석 방법: 림마, EdgeR, DESeq2
Play Video
PDF DOI DOWNLOAD MATERIALS LIST

Cite this Article

Liu, S., Wang, Z., Zhu, R., Wang,More

Liu, S., Wang, Z., Zhu, R., Wang, F., Cheng, Y., Liu, Y. Three Differential Expression Analysis Methods for RNA Sequencing: limma, EdgeR, DESeq2. J. Vis. Exp. (175), e62528, doi:10.3791/62528 (2021).

Less
Copy Citation Download Citation Reprints and Permissions
View Video

Get cutting-edge science videos from JoVE sent straight to your inbox every month.

Waiting X
Simple Hit Counter