[DADA2] plotQualityProfile()시 Error와 해결 방법

dada2 분석 중 이러한 에러 메시지를 보았다.

 

# fastq.gz이 있는 경로를 import한다 
path <- "D:/Project/DT/02.cutadapt/01" # Raw data 
list.files(path)

# 파일 경로와 파일의 이름을 합해준다 
fnFs <- sort(list.files(path, pattern="_1.fastq", full.names = TRUE))
fnRs <- sort(list.files(path, pattern="_2.fastq", full.names = TRUE))
fnFs
fnRs

# 파일 이름에서 샘플 이름을 추출한다
sample.names <- sapply(strsplit(basename(fnFs), "_1.fastq.gz"), `[`, 1)
sample.names #%>% head() 

# QC plot을 확인해 본다 
plotQualityProfile(fnFs[[1]])
# Error: BiocParallel errors
# 1 remote errors, element index: 1
# 0 unevaluated and other errors
# first remote error:
#   Error in density.default(qscore): 'x' contains missing values
# In addition: Warning message:
#   In serialize(data, node$con) :
#   'package:stats' may not be available when loading

에러메시지를 잘 살펴보자. 이는 총 두 가지를 나타낸다. 

# Error: BiocParallel errors

# In addition: Warning message: # In serialize(data, node$con) : # 'package:stats' may not be available when loading

 

 

 

일단 BiocParallel errors메시지를 검색해 본 결과 이미 논의된 두 가지의 게시글을 볼 수가 있었다.

1. https://github.com/benjjneb/dada2/issues/159

 

2. https://github.com/benjjneb/dada2/issues/164

각각 글의 해결 방법을  따라해 보자.

 

---

 

1. ShortRead 패키지의 readFastq로 원 데이터가 문제인지 확인해라 

library(ShortRead)


length(fnFs)
class(fnFs)

foo1 <- readFastq(fnFs[[1]])
foo1
quality(foo1)
foo2 <- readFastq(fnFs[[2]])
foo2
quality(foo2)

서열의 퀄리티에는 문제는 없어 보인다 

 

2. bioconductor를 업데이트해라 

if (!require("BiocManager", quietly = TRUE))
  install.packages("BiocManager")
BiocManager::install(version = "3.16")
packageVersion("BiocParallel") # 1.31.12

최신 버전은 1.32.5지만 받으려면 직접 다운로드하여야 한다.

일단 이후 plotQualityProfile()을 해봤지만 그대로 에러 메시지가 등장하였다. 

다른 방법을 찾아 보아야 한다. 

 

 

3. Cutadapt단계의 문제일 수 있다

benjjneb(dada2개발자)가 관련 글에서 이러한 말을 적었다.

"The better solution here is to enforce a minimum length using cutadapt with the --minimum-length N flag."

 

즉, cutadapt당시 최소 길이를 설정하지 않으면 길이가 0인 서열 때문에 오류가 난다는 것이었다.

구글링 해보니 적당히 20 이상 잡으면 된다고 한다. 이후 Cutadapt시 --minimum-length 옵션을 20 정도로 잡고 다시 돌렸다.

 

이후 dada2::plotQualityProfile()에서 에러가 나타나지 않았다😭

plotQualityProfile(c(fnFs[[1]], fnRs[[1]]))

 

역시 개발자의 답변이 최고다...👍

 

 

 

 

 

 

 

 

---