爾云間? 一個專門做科研的團隊

原創(chuàng)? 小果 ? 生信果?? 關(guān)注我們

今天小果繪制一下基因位置的環(huán)形圖，繪圖代碼如下：

1、安裝所需要的包

BiocManager::install(“RCircos”)

install.packages(“magritter”)

install.packages(“tidyverse”)

2、導(dǎo)入需要的R包

library(RCircos)

library(magritter)

library(tidyverse)

3、讀取數(shù)據(jù)

#圖中展示目標(biāo)基因文件

genes <- read.table("gene.txt", header = T)$SYMBOL

# 根據(jù)基因集提取出目的基因所在位置

gene_pos <- import("hg38.gtf") %>% # 載入gtf文件

? as.data.frame %>%

? # 僅選擇基因，去除轉(zhuǎn)錄本等等

? filter(source == "HAVANA", type == "gene") %>% ?

? # 保留基因位置和名稱

? dplyr::select(seqnames, start, end, gene_name) %>%? ?

? # 挑選目的基因

? filter(gene_name %in% genes)

# 原文沒說內(nèi)圈散點圖代表什么特征，我們這里隨機生成一列數(shù)值

gene_pos$gene_dot <- rnorm(nrow(gene_pos), 0, 2)

# 保存到文件，便于套用格式

write.csv(gene_pos,"gene_pos.csv", row.names = F, quote = F)





#加載基因所在的位置和數(shù)值

gene_pos <- read.csv("gene_pos.csv", header = T)

# 加載染色體Ideogram

(data("UCSC.HG38.Human.CytoBandIdeogram"))

pdf(file="circGene.pdf", height=5, width=5)

# 根據(jù)hg38構(gòu)建染色體位置，只保留chr1-22,X,Y，在圈內(nèi)部構(gòu)建三圈軌道

RCircos.Set.Core.Components(UCSC.HG38.Human.CytoBandIdeogram,

??????????????????????????? chr.exclude = NULL,

??????????????????????????? tracks.inside = 3,

??????????????????????????? tracks.outside = 0)

RCircos.Set.Plot.Area()

# 繪制染色體

RCircos.Chromosome.Ideogram.Plot()

# 在第一圈用散點在基因所在的位置標(biāo)注數(shù)值

# 調(diào)整配色

params <- RCircos.Get.Plot.Parameters()

params$track.background <- "grey" # 第三圈默認(rèn)配色為wheat，模仿原文修改為灰色

RCircos.Reset.Plot.Parameters(params)



RCircos.Scatter.Plot(gene_pos,

???????????????????? data.col = 5, # 用第5列的數(shù)值作為點的縱坐標(biāo)

???????????????????? by.fold = 1, # 點的顏色cutoff，大于等于1的基因顯示為紅色點，小于等于-1的顯示為藍(lán)色點，-1到1之間為黑點

???????????????????? track.num = 1,

???????????????????? side = "in")

# 在第二圈繪制線段標(biāo)注基因所在的位置

RCircos.Gene.Connector.Plot(genomic.data = gene_pos,

??????????????????????????? track.num = 2,

??????????????????????????? side = "in")

# 在第三圈標(biāo)注基因名

RCircos.Gene.Name.Plot(gene_pos,

?????????????????????? name.col = 4,

?????????????????????? track.num = 3,

?????????????????????? side = "in")

dev.off()

小果最終繪制出了目標(biāo)基因位置環(huán)形圖，有需要的可以借鑒學(xué)習(xí)。

推薦閱讀

• 多種花樣都能懂，富集圖的看法

• 相關(guān)性圖譜之相關(guān)性熱圖學(xué)習(xí)

• 小果教你快速優(yōu)雅使用Genecards數(shù)據(jù)庫

• 從UCSC下載TCGA數(shù)據(jù)

• 簡單高效利用Batch Entrez批量獲取基因別名

“生信果”? 生信入門、R語言、生信圖解讀與繪制、軟件操作、代碼復(fù)現(xiàn)、生信硬核知識技能、服務(wù)器、生物信息學(xué)的教程，以及基于R的分析和可視化等原創(chuàng)內(nèi)容，一起見證小白和大佬的成長。