試做你的核苷酸頻率圖
現在該更仔細檢視每個循環(cycle)的核苷酸頻率了。最好的方式是做一張視覺化圖表。通常前幾個循環會比較隨機,之後隨著循環增加,核苷酸的頻率應該會趨於穩定。
這個練習使用完整的 fastq 檔案 SRR1971253,並已替你做了一些前處理:
library(ShortRead)
fqsample <- readFastq(dirPath = "data",
pattern = "SRR1971253.fastq")
# extract reads
abc <- alphabetByCycle(sread(fqsample))
# Transpose nucleotides A, C, G, T per column
nucByCycle <- t(abc[1:4,])
# Tidy dataset
nucByCycle <- nucByCycle %>%
as_tibble() %>% # convert to tibble
mutate(cycle = 1:50) # add cycle numbers
你的任務是使用 tidyverse 函式製作一張「依循環的核苷酸頻率」折線圖!
本練習屬於課程
R 中的 Bioconductor 入門
練習說明
- 先用
glimpse()檢視nucByCycle物件的資料概況。 - 使用
pivot_longer()將核苷酸字母轉為alphabet欄位,並產生新的count欄位。 - 繪製以
cycle為 x 軸、count為 y 軸、依alphabet著色的折線圖。
動手互動練習
試著完成這個範例程式碼,體驗一下這個練習。
# Glimpse nucByCycle
___
# Create a line plot of cycle vs. count
nucByCycle %>%
# Gather the nucleotide letters in alphabet and get a new count column
pivot_longer(-cycle, names_to = ___, values_to = ___) %>%
ggplot(aes(x = ___, y = ___, color = ___)) +
geom_line(size = 0.5 ) +
labs(y = "Frequency") +
theme_bw() +
theme(panel.grid.major.x = element_blank())