• 主页
  • 课程

    关于课程

    • 课程归档
    • 成为一名讲师
    • 讲师信息
    同等学历教学

    同等学历教学

    免费
    阅读更多
  • 特色
    • 展示
    • 关于我们
    • 问答
  • 事件
  • 个性化
  • 博客
  • 联系
  • 站点资源
    有任何问题吗?
    (00) 123 456 789
    weinfoadmin@weinformatics.cn
    注册登录
    恒诺新知
    • 主页
    • 课程

      关于课程

      • 课程归档
      • 成为一名讲师
      • 讲师信息
      同等学历教学

      同等学历教学

      免费
      阅读更多
    • 特色
      • 展示
      • 关于我们
      • 问答
    • 事件
    • 个性化
    • 博客
    • 联系
    • 站点资源

      未分类

      • 首页
      • 博客
      • 未分类
      • circlize 之可视化基因组数据

      circlize 之可视化基因组数据

      • 发布者 weinfoadmin
      • 分类 未分类, 老俊俊的生信笔记
      • 日期 2021年9月10日
      • 评论 0评论

      感谢老俊俊的大力支持。我们会每日跟新,欢迎您关注老俊俊的生信笔记。


      点击上方公众号查看更多精彩内容



      circlize 可用来在基因组学研究领域可视化及展示基因组之间不同染色体的关联性。

      为了便于基因组学分析,circlize 特别提供了专注于绘制基因组图的功能。这些函数与基本图形函数是类似的,但是需要输入特殊格式的数据:

      参数 说明
      circos.genomicTrack() 创建轨道和添加图层
      circos.genomicPoints() 添加点
      circos.genomicLines() 添加线
      circos.genomicRect() 添加矩形
      circos.genomicText() 添加文字
      circos.genomicLink() 添加连接

      输入数据:

      基因组数据通常存储为一个 table 形式,其中前三列定义基因组区域,随后的列是与相应区域相关联的值。每个基因组区域由三个要素组成:基因组范畴(多数为染色体)、基因组范畴上的起始位置和结束位置。这种数据结构被称为 BED 格式,广泛应用于基因组研究。

      circlize 提供了一个简单的函数 generateRandomBed(),它生成随机的基因组数据。位置是由人类基因组均匀地产生的,染色体上的区域数与染色体的长度近似成正比。其中,nr 和 nc 参数控制用户需要的行数和列数。请注意,nr 与函数返回的行数不完全相同。Fun 参数是用来生成随机值的自定义函数。

      set.seed(999)
      bed = generateRandomBed()
      head(bed)
      ##    chr   start     end       value1
      ## 1 chr1  118750  451929 -0.367702502
      ## 2 chr1  472114  805024 -0.001764946
      ## 3 chr1  807013  914103 -0.668379255
      ## 4 chr1 1058915 1081590  0.390291705
      ## 5 chr1 1194426 1341960 -1.305901261
      ## 6 chr1 1670402 2048723  0.340227447

      生成特定行数和列数:

      bed = generateRandomBed(nr = 200, nc = 4)
      nrow(bed)
      ## [1] 205

      自定义函数添加列:

      bed = generateRandomBed(nc = 2, fun = function(k) sample(letters, k, replace = TRUE))
      head(bed)
      ##    chr   start     end value1 value2
      ## 1 chr1  275067  349212      q      e
      ## 2 chr1  350892  658970      u      f
      ## 3 chr1  674620  875538      y      r
      ## 4 chr1 1053255 1115056      p      e
      ## 5 chr1 1127139 1545066      a      s
      ## 6 chr1 1864092 1973553      x      u


      1、cytoband 基因组数据初始化


      cytoband 数据是初始化基因组图的理想数据源。它包含染色体的长度以及所谓的染色体带注释,以帮助识别染色体上的位置。

      使用 circos.initializeWithIdeogram()函数加载 cytoband 基因组数据,默认加载的是 hg19 的基因组数据:

      circos.initializeWithIdeogram()
      text(0, 0, "default", cex = 1)

      查看信息:

      circos.info()
      ## All your sectors:
      ##  [1] "chr1"  "chr2"  "chr3"  "chr4"  "chr5"  "chr6"  "chr7"  "chr8"  "chr9"
      ## [10] "chr10" "chr11" "chr12" "chr13" "chr14" "chr15" "chr16" "chr17" "chr18"
      ## [19] "chr19" "chr20" "chr21" "chr22" "chrX"  "chrY"
      ##
      ## All your tracks:
      ## [1] 1 2
      ##
      ## Your current sector.index is chrY
      ## Your current track.index is 2
      circos.clear()

      也可以通过指定物种参数初始化其他物种,并自动下载对应物种的 cytoband 数据文件:

      circos.initializeWithIdeogram(species = "hg18")
      circos.initializeWithIdeogram(species = "mm10")

      当你在处理稀有物种和没有可用的 cytoband 数据时,circos.initializeWithIdeogram()将尝试继续下载 UCSC 的 chromInfo 文件,其中也包含染色体的长度,但是当然,在图上没有 ideogram track。

      在某些情况下,当没有互联网连接或 UCSC 上没有相应的数据可用的。可以手动构造一个包含染色体范围的 table,如果它存储在一个文件中或目录里,并发送到 circos.initializeWithIdeogram():

      # 放在目录下
      cytoband.file = system.file(package = "circlize", "extdata", "cytoBand.txt")
      circos.initializeWithIdeogram(cytoband.file)

      # 读入数据
      cytoband.df = read.table(cytoband.file, colClasses = c("character", "numeric",
          "numeric", "character", "character"), sep = "t")
      circos.initializeWithIdeogram(cytoband.df)

      默认情况下,circos.initiizewithideogram()使用数据中可用的所有染色体来初始化环形图。用户可以通过指定 chromosome.index 来选择染色体的一个子集。这个参数也适用于染色体排序:

      circos.initializeWithIdeogram(chromosome.index = paste0("chr", c(3,5,2,8)))
      text(0, 0, "subset of chromosomes", cex = 1)
      circos.clear()

      在初始化环形图之后,circos.initializeWithIdeogram()另外创建一个轨道,其中有基因组轴和染色体名称,并在有 ideogram 的地方创建另一个轨道(取决于 cytoband 数据是否可用)。plotType 参数用于控制添加哪种类型的轨道:

      circos.initializeWithIdeogram(plotType = c("axis", "labels"))
      text(0, 0, "plotType = c('axis', 'labels')", cex = 1)
      circos.clear()

      circos.initializeWithIdeogram(plotType = NULL)
      text(0, 0, "plotType = NULL", cex = 1)
      circos.clear()

      我们可以用 circos.par()函数来设置其它参数,需要 circos.clear()来完成绘图:

      circos.par("start.degree" = 90)
      circos.initializeWithIdeogram()
      circos.clear()
      text(0, 0, "'start.degree' = 90", cex = 1)

      circos.par("gap.degree" = rep(c(2, 4), 12))
      circos.initializeWithIdeogram()
      circos.clear()
      text(0, 0, "'gap.degree' = rep(c(2, 4), 12)", cex = 1)


      2、自定义染色体轨道


      默认情况下,circos.initializeWithIdeogram()初始化布局并添加两条轨道。当 plotType 参数设置为 NULL 时,只初始化环形布局,但不添加任何内容。这使得用户可以完全设计自己风格的染色体轨道。

      在下面的例子中,我们使用不同的颜色来代表染色体,并将染色体的名称放在每个单元格的中心:

      set.seed(123)
      circos.initializeWithIdeogram(plotType = NULL)
      circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
          chr = CELL_META$sector.index
          xlim = CELL_META$xlim
          ylim = CELL_META$ylim
          circos.rect(xlim[1], 0, xlim[2], 1, col = rand_color(1))
          circos.text(mean(xlim), mean(ylim), chr, cex = 0.7, col = "white",
              facing = "inside", niceFacing = TRUE)
      }, track.height = 0.15, bg.border = NA)
      circos.clear()



      3、初始化其它基因组类别


      染色体只是基因组范畴的一种特殊情况。circos.genomicInitialize()可以用任何类型的基因组类别初始化环形布局。事实上,circos.initializeWithIdeogram()是由 circos.genomicInitialize()实现的。circos.genomicInitialize()的输入数据也是至少有三列的数据。第一列是基因组类别(对于细胞带数据,是染色体名称),接下来的两列是每个基因组类别中的位置。每个类别的范围将被推断为对应类别的最小位置和最大位置。

      在下面的例子中,是一个由三个基因初始化环形图的,构造一个数据:

      df = data.frame(
          name  = c("TP53",  "TP63",    "TP73"),
          start = c(7565097, 189349205, 3569084),
          end   = c(7590856, 189615068, 3652765))
      circos.genomicInitialize(df)

      在以下示例中,我们在圆形布局中绘制 TP53,TP63 和 TP73 的转录本:

      tp_family = readRDS(system.file(package = "circlize", "extdata", "tp_family_df.rds"))
      head(tp_family)
      ##   gene   start     end        transcript exon
      ## 1 TP53 7565097 7565332 ENST00000413465.2    7
      ## 2 TP53 7577499 7577608 ENST00000413465.2    6
      ## 3 TP53 7578177 7578289 ENST00000413465.2    5
      ## 4 TP53 7578371 7578554 ENST00000413465.2    4
      ## 5 TP53 7579312 7579590 ENST00000413465.2    3
      ## 6 TP53 7579700 7579721 ENST00000413465.2    2

      创建 3 个轨道:

      circos.genomicInitialize(tp_family)
      circos.track(ylim = c(0, 1),
          bg.col = c("#FF000040", "#00FF0040", "#0000FF40"),
          bg.border = NA, track.height = 0.05)

      给每个转录本绘制线和矩形连接起来:

      n = max(tapply(tp_family$transcript, tp_family$gene, function(x) length(unique(x))))
      circos.genomicTrack(tp_family, ylim = c(0.5, n + 0.5),
          panel.fun = function(region, value, ...) {
              all_tx = unique(value$transcript)
              for(i in seq_along(all_tx)) {
                  l = value$transcript == all_tx[i]
                  # for each transcript
                  current_tx_start = min(region[l, 1])
                  current_tx_end = max(region[l, 2])
                  circos.lines(c(current_tx_start, current_tx_end),
                      c(n - i + 1, n - i + 1), col = "#CCCCCC")
                  circos.genomicRect(region[l, , drop = FALSE], ytop = n - i + 1 + 0.4,
                      ybottom = n - i + 1 - 0.4, col = "orange", border = NA)
              }
      }, bg.border = NA, track.height = 0.4)
      circos.clear()

      你可能会注意到坐标轴的开始变成了 0KB,而不是原来的值。这只是一个轴标签的调整,以反映到每个基因开始的相对距离,而细胞中的坐标仍然使用原始值。设置 tickLabelsStartFromZero 为 FALSE 将坐标轴标签恢复为原始值。


      4、缩放染色体


      我们首先定义一个函数 extend.chromosome(),它将染色体子集中的数据复制到原始数据中:

      extend_chromosomes = function(bed, chromosome, prefix = "zoom_") {
          zoom_bed = bed[bed[[1]] %in% chromosome, , drop = FALSE]
          zoom_bed[[1]] = paste0(prefix, zoom_bed[[1]])
          rbind(bed, zoom_bed)
      }

      我们使用 read.cytoband()从 UCSC 下载并读取 cytoband 数据。接下来,正常染色体和放大染色体的 x 范围分别归一化:

      cytoband = read.cytoband()
      cytoband_df = cytoband$df
      chromosome = cytoband$chromosome

      xrange = c(cytoband$chr.len, cytoband$chr.len[c("chr1", "chr2")])
      normal_chr_index = 1:24
      zoomed_chr_index = 25:26

      # normalize in normal chromsomes and zoomed chromosomes separately
      sector.width = c(xrange[normal_chr_index] / sum(xrange[normal_chr_index]),
                       xrange[zoomed_chr_index] / sum(xrange[zoomed_chr_index]))

      绘图:

      circos.par(start.degree = 90)
      circos.initializeWithIdeogram(extend_chromosomes(cytoband_df, c("chr1", "chr2")),
          sector.width = sector.width)

      添加新的轨道:

      bed = generateRandomBed(500)
      circos.genomicTrack(extend_chromosomes(bed, c("chr1", "chr2")),
          panel.fun = function(region, value, ...) {
              circos.genomicPoints(region, value, pch = 16, cex = 0.3)
      })

      添加连接:

      circos.link("chr1", get.cell.meta.data("cell.xlim", sector.index = "chr1"),
          "zoom_chr1", get.cell.meta.data("cell.xlim", sector.index = "zoom_chr1"),
          col = "#00000020", border = NA)
      circos.clear()



      5、整合两个基因组


      在某些情况下,希望在环形图中可视化多个基因组。这可以通过合成一个组合基因组来实现。在下面的例子中,结合了人类和老鼠的基因组:

      human_cytoband = read.cytoband(species = "hg19")$df
      mouse_cytoband = read.cytoband(species = "mm10")$df

      # 添加前缀来区分
      human_cytoband[ ,1] = paste0("human_", human_cytoband[, 1])
      mouse_cytoband[ ,1] = paste0("mouse_", mouse_cytoband[, 1])

      合并:

      cytoband = rbind(human_cytoband, mouse_cytoband)
      head(cytoband)
      ##           V1       V2       V3     V4     V5
      ## 1 human_chr1        0  2300000 p36.33   gneg
      ## 2 human_chr1  2300000  5400000 p36.32 gpos25
      ## 3 human_chr1  5400000  7200000 p36.31   gneg
      ## 4 human_chr1  7200000  9200000 p36.23 gpos25
      ## 5 human_chr1  9200000 12700000 p36.22   gneg
      ## 6 human_chr1 12700000 16200000 p36.21 gpos50

      排序,使让人类 1 号染色体接近老鼠 1 号染色体:

      chromosome.index = c(paste0("human_chr", c(1:22, "X", "Y")),
                           rev(paste0("mouse_chr", c(1:19, "X", "Y"))))
      circos.initializeWithIdeogram(cytoband, chromosome.index = chromosome.index)
      circos.clear()

      关掉染色体名称和轴,新建一个小轨道来区分人类染色体和老鼠染色体用 highlight.chromosome()函数,对于染色体,我只写数值索引(也包括 X 和 Y),人类和老鼠的染色体间隔 5。(circos.par("gap.after")):

      circos.par(gap.after = c(rep(1, 23), 5, rep(1, 20), 5))
      circos.initializeWithIdeogram(cytoband, plotType = NULL,
          chromosome.index = chromosome.index)
      circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
          circos.text(CELL_META$xcenter, CELL_META$ylim[2] + mm_y(2),
              gsub(".*chr", "", CELL_META$sector.index), cex = 0.6, niceFacing = TRUE)
      }, track.height = mm_h(1), cell.padding = c(0, 0, 0, 0), bg.border = NA)
      highlight.chromosome(paste0("human_chr", c(1:22, "X", "Y")),
          col = "red", track.index = 1)
      highlight.chromosome(paste0("mouse_chr", c(1:19, "X", "Y")),
          col = "blue", track.index = 1)

      circos.genomicIdeogram(cytoband)
      circos.clear()

      也只能通过染色体范围,即每个染色体的长度来创建。在下面的代码中,read.chromInfo()可以获取特定基因组的染色体范围:

      human_chromInfo = read.chromInfo(species = "hg19")$df
      mouse_chromInfo = read.chromInfo(species = "mm10")$df
      human_chromInfo[ ,1] = paste0("human_", human_chromInfo[, 1])
      mouse_chromInfo[ ,1] = paste0("mouse_", mouse_chromInfo[, 1])
      chromInfo = rbind(human_chromInfo, mouse_chromInfo)
      # note the levels of the factor controls the chromosome orders in the plot
      chromInfo[, 1] = factor(chromInfo[ ,1], levels = chromosome.index)
      head(chromInfo)
      ##          chr start       end
      ## 1 human_chr1     0 249250621
      ## 2 human_chr2     0 243199373
      ## 3 human_chr3     0 198022430
      ## 4 human_chr4     0 191154276
      ## 5 human_chr5     0 180915260
      ## 6 human_chr6     0 171115067

      对于指定的染色体:范围,使用 circos.genomicInitialize()来初始化布局:

      circos.par(gap.after = c(rep(1, 23), 5, rep(1, 20), 5))
      circos.genomicInitialize(chromInfo, plotType = NULL)
      circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
          circos.text(CELL_META$xcenter, CELL_META$ylim[2] + mm_y(2),
              gsub(".*chr", "", CELL_META$sector.index), cex = 0.6, niceFacing = TRUE)
      }, track.height = mm_h(1), cell.padding = c(0, 0, 0, 0), bg.border = NA)
      highlight.chromosome(paste0("human_chr", c(1:22, "X", "Y")),
          col = "red", track.index = 1)
      highlight.chromosome(paste0("mouse_chr", c(1:19, "X", "Y")),
          col = "blue", track.index = 1)

      circos.track(ylim = c(0, 1))
      circos.clear()

      给人和小鼠添加点和连接:

      circos.par(gap.after = c(rep(1, 23), 5, rep(1, 20), 5))
      circos.genomicInitialize(chromInfo, plotType = NULL)
      circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
          circos.text(CELL_META$xcenter, CELL_META$ylim[2] + mm_y(2),
              gsub(".*chr", "", CELL_META$sector.index), cex = 0.6, niceFacing = TRUE)
      }, track.height = mm_h(1), cell.padding = c(0, 0, 0, 0), bg.border = NA)
      highlight.chromosome(paste0("human_chr", c(1:22, "X", "Y")),
          col = "red", track.index = 1)
      highlight.chromosome(paste0("mouse_chr", c(1:19, "X", "Y")),
          col = "blue", track.index = 1)

      circos.genomicIdeogram(cytoband)

      # a track of points
      human_df = generateRandomBed(200, species = "hg19")
      mouse_df = generateRandomBed(200, species = "mm10")
      human_df[ ,1] = paste0("human_", human_df[, 1])
      mouse_df[ ,1] = paste0("mouse_", mouse_df[, 1])
      df = rbind(human_df, mouse_df)
      circos.genomicTrack(df, panel.fun = function(region, value, ...) {
          circos.genomicPoints(region, value, col = rand_color(1), cex = 0.5, ...)
      })

      # links between human and mouse genomes
      human_mid = data.frame(
          chr = paste0("human_chr", 1:19),
          mid = round((human_chromInfo[1:19, 2] + human_chromInfo[1:19, 3])/2)
      )
      mouse_mid = data.frame(
          chr = paste0("mouse_chr", 1:19),
          mid = round((mouse_chromInfo[1:19, 2] + mouse_chromInfo[1:19, 3])/2)
      )
      circos.genomicLink(human_mid, mouse_mid, col = rand_color(19))
      circos.clear()
      text(-0.9, -0.8, "Humanngenome")
      text(0.9, 0.8, "Mousengenome")



      end


      发现更多精彩

      关注公众号

      欢迎小伙伴留言评论!

      点击我留言!

      今天的分享就到这里了,敬请期待下一篇!

      最后欢迎大家分享转发,您的点赞是对我的鼓励和肯定!

      如果觉得对您帮助很大,打赏一下吧!

      请关注“恒诺新知”微信公众号,感谢“R语言“,”数据那些事儿“,”老俊俊的生信笔记“,”冷🈚️思“,“珞珈R”,“生信星球”的支持!

      • 分享:
      作者头像
      weinfoadmin

      上一篇文章

      circlize 之 Create plotting regions
      2021年9月10日

      下一篇文章

      circlize 之 Advanced layout
      2021年9月10日

      你可能也喜欢

      2-1675088548
      lncRNA和miRNA生信分析系列讲座免费视频课和课件资源包,干货满满
      30 1月, 2023
      9-1675131201
      如何快速批量修改 Git 提交记录中的用户信息
      26 1月, 2023
      8-1678501786
      肿瘤细胞通过改变CD8+ T细胞中的丙酮酸利用和琥珀酸信号来调控抗肿瘤免疫应答。
      7 12月, 2022

      留言 取消回复

      要发表评论,您必须先登录。

      搜索

      分类

      • R语言
      • TCGA数据挖掘
      • 单细胞RNA-seq测序
      • 在线会议直播预告与回放
      • 数据分析那些事儿分类
      • 未分类
      • 生信星球
      • 老俊俊的生信笔记

      投稿培训

      免费

      alphafold2培训

      免费

      群晖配置培训

      免费

      最新博文

      Nature | 单细胞技术揭示衰老细胞与肌肉再生
      301月2023
      lncRNA和miRNA生信分析系列讲座免费视频课和课件资源包,干货满满
      301月2023
      如何快速批量修改 Git 提交记录中的用户信息
      261月2023
      logo-eduma-the-best-lms-wordpress-theme

      (00) 123 456 789

      weinfoadmin@weinformatics.cn

      恒诺新知

      • 关于我们
      • 博客
      • 联系
      • 成为一名讲师

      链接

      • 课程
      • 事件
      • 展示
      • 问答

      支持

      • 文档
      • 论坛
      • 语言包
      • 发行状态

      推荐

      • iHub汉语代码托管
      • iLAB耗材管理
      • WooCommerce
      • 丁香园论坛

      weinformatics 即 恒诺新知。ICP备案号:粤ICP备19129767号

      • 关于我们
      • 博客
      • 联系
      • 成为一名讲师

      要成为一名讲师吗?

      加入数以千计的演讲者获得100%课时费!

      现在开始

      用你的站点账户登录

      忘记密码?

      还不是会员? 现在注册

      注册新帐户

      已经拥有注册账户? 现在登录

      close
      会员购买 你还没有登录,请先登录
      • ¥99 VIP-1个月
      • ¥199 VIP-半年
      • ¥299 VIP-1年
      在线支付 激活码

      立即支付
      支付宝
      微信支付
      请使用 支付宝 或 微信 扫码支付
      登录
      注册|忘记密码?