背景介绍

细胞互作图

这张图来自于一篇对胎盘母胎界面的细胞互作研究[1]，这篇文献筛选出了所有细胞表达的配体和受体，利用现有的数据库找到配体-受体对，用箭头将这些细胞表达对应的配体-受体对连接起来，从而推断出不同类型细胞间的互作关系。
从外至内的第二圈7个色块展示的是细胞类型，最外圈展示它们所表达的配体和受体，红绿色块表示该基因的表达水平，黑色箭头表示母胎细胞间的互作，灰色箭头指的是细胞滋养层间的互作。

如果你对这种图感兴趣的话，那就跟我一起画一画吧！

目的

利用细胞表达的配体和受体，将它们的互作关系模型用圆圈图展示出来。

方法

1. 读取文件

df = read.csv('cordblood.csv',header = T)

这个文件共有5列，三个分类标准：分别是第一列细胞类型cell；第二列配体还是受体lr；第三列基因名gene。最后一列fc为我们数据的具体数值foldchange值，也就是差异倍数。gene_id只起一个序号的作用，方便之后画图使用。

cordblood.csv 示前25行

2. 整体布局与初始化

首先加载包，调整画布范围，避免基因标签出界（默认画布大小为c(1,1)c(1,1)）。circlize的画图逻辑一定是从外至内一圈圈地叠加，并且最外圈的半径一定是1，因此若要在最外圈的外面添加内容的话，我们的画布半径要比1大，这里设置为1.1。

library(circlize) 
# circos.clear()  #这个命令用于清空画布，画错时要运行此命令重新再画。
#整体布局
circos.par(canvas.xlim =c(-1.1,1.1),canvas.ylim = c(-1.1,1.1),cell.padding = c(0.02,0,0.02,0))

3. 画第一圈

画第一圈一定要使用circos.initialize进行初始化，初始化需要设定factors。factors可以粗暴地理解为最外圈有多少个方块种类。我们这里要将每个基因的fc值用色块表示出来，有多少行就画多少个小色块，因此使用的factor就是gene_id。这里为什么不用gene呢？是因为不同的细胞中会表达相同类型的基因，如此一来得到的factor就少了。

fa = df$gene_id
fa = factor(fa,levels = fa)
circos.initialize(factors = fa, xlim = c(0,1)) # 初始化

这一行运行完之后你还不会看到任何东西，因为它只是在画布上设定了一个轨道（track），把每一块位置计算分割好了，并没有填充内容，下面我们来将这一圈填充内容。
circos.trackPlotRegion是绘画的关键，trackPlotRegion顾名思义就是在轨道上的区域上作画。参数说明：
track.height设置圆圈环的宽度，比如说0.15就表示这是一个外圈半径为1内圈半径为0.85的环（前面说了最外圈半径一定为1）。
bg.border设置每个区域的边界颜色。
bg.col设置每个区域的填充颜色，这里填充的是根据fc值大小所对应的 黑-黄-红 过渡颜色。

# 设置fc的大小对应的颜色，随着颜色从黑到黄到红过渡，fc值从-10至0至10
col_fun = colorRamp2(c(-10, 0, 10), c("black", "yellow", "red"))

circos.trackPlotRegion(
  ylim = c(0, 1), track.height = 0.15, bg.border = 'black', bg.col = col_fun(df$fc),
  panel.fun = function(x, y) {
    sector.index = get.cell.meta.data('sector.index')
    xlim = get.cell.meta.data('xlim')
    ylim = get.cell.meta.data('ylim')
  } )

给每个区域加上标签——它们的基因名。用circos.axis命令设置label。
sector.index即我们的区域索引第四列gene_id，
labels指的是具体要标注上去的内容，也就是第三列基因名。
其他的参数是调整标签的字体大小（labels.cex）、颜色（col）、方向（labels.facing）之类的。

# 标注基因
for(i in 1:nrow(df)){
  circos.axis(sector.index= df[i,4], direction = "outside", labels=df[i,3], 
              labels.facing = "clockwise",labels.cex=.58, col = 'black',
              labels.away.percentage=0.1, minor.ticks=0, major.at=seq(1, length(df$gene)))
}

4.画第二圈 细胞类型

还是用circos.trackPlotRegion进行绘画，圆圈环的宽度还是0.15。根据每个细胞类型对应的gene_id号进行上色，用highlight.sector给它们进行高亮以及标注。这里注意 track.index = 2，选择高亮第二圈。font设置字体大小，col为区域填充色。niceFacing为标注字体的方向是不是方便人类阅读。大家可以自行尝试一下niceFacing = T会有什么不同的样子。

# 第二圈 细胞类型
circos.trackPlotRegion(
  ylim = c(0, 1), track.height = 0.15, bg.border = NA,
  panel.fun = function(x, y) {
    sector.index = get.cell.meta.data('sector.index')
    xlim = get.cell.meta.data('xlim')
    ylim = get.cell.meta.data('ylim')
  } )

# cell
highlight.sector(as.character(df$gene_id[1:11]), track.index = 2,
                 text = 'B cell', niceFacing = F, font = 2, col = '#CCEBC5')
highlight.sector(as.character(df$gene_id[12:75]), track.index = 2,
                 text = 'mon,neu,eos ', niceFacing = F, font = 2, col ='#FFFFB3')
highlight.sector(as.character(df$gene_id[76:142]), track.index = 2,
                 text = 'others', niceFacing = F, font = 2)
highlight.sector(as.character(df$gene_id[143:170]), track.index = 2,
                 text = 'T cell', niceFacing = F, font = 2, col = '#FFCC99')
highlight.sector(as.character(df$gene_id[171:196]), track.index = 2,
                 text = 'NK', niceFacing = F, font = 2, col = '#ccccff')

5. 画第三圈 配体与受体

这次设置的track.height为0.08，窄一点。这一圈的绘画原理同第二圈，找出每个区域对应的编号范围，进行高亮。text.col设置标注颜色。代码看上去很长，但其实找到规律看一行就够了。

#第三圈 配体受体
circos.trackPlotRegion(
  ylim = c(0, 1), track.height = 0.08, bg.border = NA,
  panel.fun = function(x, y) {
    sector.index = get.cell.meta.data('sector.index')
    xlim = get.cell.meta.data('xlim')
    ylim = get.cell.meta.data('ylim')
  } )
#红蓝配
highlight.sector(as.character(df$gene_id[1:6]), track.index = 3,
                 text = 'ligand', niceFacing = F,col = '#FB8072',text.col = 'white')
highlight.sector(as.character(df$gene_id[7:11]), track.index = 3,
                 text = 'receptor', niceFacing = F,col = '#80B1D3',text.col = 'white')
highlight.sector(as.character(df$gene_id[12:43]), track.index = 3,
                 text = 'ligand', niceFacing = F,col = '#FB8072',text.col = 'white')
highlight.sector(as.character(df$gene_id[44:75]), track.index = 3,
                 text = 'receptor', niceFacing = F,col = '#80B1D3',text.col = 'white')
highlight.sector(as.character(df$gene_id[76:110]), track.index = 3,
                 text = 'ligand', niceFacing = F,col = '#FB8072',text.col = 'white')
highlight.sector(as.character(df$gene_id[111:142]), track.index = 3,
                 text = 'receptor', niceFacing = F,col = '#80B1D3',text.col = 'white')
highlight.sector(as.character(df$gene_id[143:155]), track.index = 3,
                 text = 'ligand', niceFacing = F,col = '#FB8072',text.col = 'white')
highlight.sector(as.character(df$gene_id[156:170]), track.index = 3,
                 text = 'receptor', niceFacing = F,col = '#80B1D3',text.col = 'white')
highlight.sector(as.character(df$gene_id[171:178]), track.index = 3,
                 text = 'ligand', niceFacing = F,col = '#FB8072',text.col = 'white')
highlight.sector(as.character(df$gene_id[179:196]), track.index = 3,
                 text = 'receptor', niceFacing = F,col = '#80B1D3',text.col = 'white')

6. 画互作关系-箭头

来到了最麻烦的一步了！但是胜利就在前方！
首先要找出这些基因中的配体-受体对的关系，才能用箭头将它们连起来。步骤就是：
（1）将数据库中的配体-受体对关系文件下载下来
（2）mapping到我们自己的数据中，找出基因中有相互作用的配体-受体对
（3）画箭头

假设你已经下载好了，现在开始处理文件：

将数据库的配体-受体对mapping到自己数据中

下面这一段代码都是处理文件的过程，用的都是非常基础的R语言数据处理知识，看不懂的童鞋自己去问谷歌吧。

## 读取数据库中配体-受体对关系文件
lr = read.csv('./ligand_receptors.csv',header = TRUE)
## 将配体和受体信息分开
ldf = subset(df,lr == 'ligand')
rdf = subset(df,lr == 'receptor')
## 创建新的数据框
lrid = lr 
lrid = within(lrid,{geneid1 = NA})
lrid = within(lrid,{geneid2 = NA})
lrid = lrid[,c('Ligand','geneid1','Receptor','geneid2')]
## 添加geneid，用于后面画箭头的索引使用。
for(i in 1:nrow(rdf)){
  index = which(lrid$Receptor %in% rdf[i,3])
  for(ind in index){lrid[ind,4] = as.character(rdf[i,4])}
}

for(i in 1:nrow(ldf)){
  index = which(lrid$Ligand %in% ldf[i,3])
  for(ind in index){lrid[ind,2] = as.character(ldf[i,4])}
}
#去除NA
lrid = na.omit(lrid)

以上一番操作后，得到的lrid数据框长这样：

处理好数据之后，我们来画箭头。利用circos.link进行绘图，前四个参数是必不可少的：
（1）sector.index1是指箭头起点的索引。
（2）point1是指起点的宽度，因为我们是一条线，所以宽度为0。如果设置为2，那么这条线就会变成一条带。
（3）sector.index2 = lrid[i,6]是指箭头终点的索引。
（4）point2是指终点的宽度，同理，线是没有宽度的。
我们要让箭头从配体指向受体，因此索引起点就是geneid1，终点就是geneid2.
（5）lty是指连接样式，1为实线，2为虚线。

#画图
for(i in 1:nrow(lrid)){
  circos.link(sector.index1 = lrid[i,3], point1 = 0, sector.index2 = lrid[i,6], point2 = 0 ,directional = 1,
              h=0.5, lwd=1, col="black",lty=1) }

当当当当！大功告成！

福利：

学习光看是不行的，真正要化为自己的知识当然得实操一下。因此我将这篇文章的代码以及所有输入文件已上传至【百度云盘】，方便大家学习使用：
链接: https://pan.baidu.com/s/1APnw2K3N05JtgV3GcBdFXg 提取码: 59yr

若要系统学习circlize包的使用方法，更多内容移步：https://jokergoo.github.io/circlize_book/book/

参考文献：

[1] Pavličev, Mihaela, et al. "Single-cell transcriptomics of the human placenta: inferring the cell communication network of the maternal-fetal interface." Genome research 27.3 (2017): 349-361.

---

文章已发布到微信公众号：百味科研芝士，欢迎关注。