#Create DESeqDataSet for DESeq2
Bulk seq analysis는 보통 Deseq2 라는 함수를 사용하여 수행합니다. 이는 DeseqDataSet(dds)라고 하는 data frame 형태의 데이터를 사용하기 때문에 이에 맞춰 데이터를 따로 만들어줘야 합니다.
Step 1. Create tx2gene for tximport
tximport라는 함수는 전사체 수준의 발현량 데이터를 유전자 수준의 발현량 데이터로 변환하는 데 사용하는 함수입니다.
이는 tx2gene이라고 하는 Transcript id 와 Gene symbol 데이터를 포함하는 dataframe 구조의 데이터를 사용합니다.
그래서 먼저 tx2gene 데이터를 만들어주는 것이 필요합니다.
우선, 참조 유전체 서열 데이터에 저장된 Transcript id와 Gene symbol만 추출하여서 TXNAME과 SYMBOL 데이터를 만들어줍니다.
#Creat TXNAME data
(Kallisto) daehwankim:~/00_Reference$ grep ">" ./Referce.Transcript.Sequences.fa > TXNAME
(Kallisto) daehwankim:~/00_Reference$ awk '{sub(/^./,"")}1' TXNAME > TXNAME.tmp && mv TXNAME.tmp TXNAME
(Kallisto) daehwankim:~/00_Reference$ awk -F "|" '{print $1}' TXNAME > TXNAME.tmp && mv TXNAME.tmp TXNAME
#Creat SYMBOL data
(Kallisto) daehwankim:~/00_Reference$ grep ">" ./Referce.Transcript.Sequences.fa > SYMBOL
(Kallisto) daehwankim:~/00_Reference$ awk '{sub(/^./,"")}1' SYMBOL > SYMBOL.tmp && mv SYMBOL.tmp SYMBOL
(Kallisto) daehwankim:~/00_Reference$ awk -F "|" '{print $6}' SYMBOL > SYMBOL.tmp && mv SYMBOL.tmp SYMBOL그 다음, 터미널을 활용하여 서버에서 로컬로 TXNAME과 SYMBOL 데이터를 export 합니다.
이때도 이전과 마찬가지로 scp 코드를 사용하면 됩니다.
(base) daehwankim-MacBookAir ~ % scp daehwankim@100.000.000.000:/home/daehwankim/00_Reference/TXNAME /Users_Path_to/Bulk_seq_practice마지막으로 R studio 상에서 read.delim 함수를 사용하여 로컬에 있는 데이터를 import 하고, 이를 dataframe 형태로 변환합니다.
library(dplyr)
TXNAME <- read.delim('/Users_Path_to/Bulk_seq_practice/TXNAME', header=F)
SYMBOL <- read.delim('/Users_Path_to/Bulk_seq_practice/SYMBOL', header=F)
tx2gene <- data.frame(TXNAME, SYMBOL)
tx2gene <- rename(tx2gene, TXNAME = V1, SYMBOL = V1.1)
Step 2. Create file path for tximport
tximport 함수는 tx2gene 외에도 kallisto를 이용하여 만든 전처리된 데이터들도 필요로 합니다.
때문에 로컬에 저장되어 있는 데이터들을 R studio에서 tximport를 통해 이용할 수 있도록 데이터들이 저장되어 있는 파일경로를 filepath 함수를 사용하여 지정하여 줍니다.
sample1 <- c('cont1.tsv', 'cont2.tsv', 'cont3.tsv','Cachexia1.tsv','Cachexia2.tsv','Cachexia3.tsv')
sample1 <- as.data.frame(sample1)
files1 <- file.path('/Users_Path_to/Cont', sample1$sample1[1:3])
files1 <- c(files1, file.path('/Users_Path_to/Cachexia', sample1$sample1[4:6]))
names(files1) <- c('Cont1.tsv' , 'Cont2.tsv' , 'Cont3.tsv' , 'Cachexia1.tsv' , 'Cachexia2.tsv' , 'Cachexia3.tsv')파일경로가 제대로 지정이 되었다면, file.exists()를 사용하여 확인했을 때 결과가 TRUE로 나오게 됩니다.
file.exists(files1)
#결과값
> file.exists(files1)
[1] TRUE TRUE TRUE TRUE TRUE TRUEStep 3. Transformation of transcript-level data into gene-level data using tximport
tx2gene과 각 tsv file의 경로를 내포하는 데이터(files1)가 모두 준비되었다면, tximport 함수를 사용하여 유전자 수준의 데이터로 변환합니다.
library(tximport)
Cont_Cachexia.tximport.tsv <- tximport(files1, type = 'kallisto', tx2gene = tx2gene, ignoreAfterBar = TRUE)
head(Cont_Cachexia.tximport.tsv$counts)Step 4. Create DESeqDataSet (dds)
tximport로 변환된 데이터를 DESeqDataSetFromTximport를 통해 다시 DESeq2 함수에 최적화된 데이터로 변환합니다.
여기서 sampleTable 데이터는 metadata로써 '데이터에 대한 정보'를 담은 데이터입니다.
즉, 각 샘플 혹은 실험에 대한 '추가적인 정보'를 담은 데이터라고 생각하시면 됩니다.
저는 이 metadata를 이용하여 Control과 Cachexia라고 하는 두가지 실험 조건이 있다는 정보를 dds에 넣어주었습니다.
library(DESeq2)
sampleTable <- data.frame(condition=factor(rep(c("Cont","Cachexia"), each = 3)))
rownames(sampleTable) <- colnames(Cont_Cachexia.tximport.tsv$counts)
dds <- DESeqDataSetFromTximport(Cont_Cachexia.tximport.tsv, sampleTable, ~condition)이번 글에서는 DEG 분석을 위해 필요한 데이터를 구축해보았습니다.
다음 글에서는 본격적으로 DESeq2 함수를 이용하여 DEG 분석을 진행해보도록 하겠습니다
감사합니다.
