在這篇短文中，您可以學(xué)習(xí)如何在Seurat v3對(duì)象上執(zhí)行基本的NicheNet分析，以及如何在circos圖中可視化輸出白指。值得注意的是，我們作為NicheNet的開發(fā)者疚宇，通常推薦通過(guò)結(jié)合幾個(gè)熱圖(配體活性、配體-靶標(biāo)鏈接赏殃、配體-受體鏈接敷待、配體表達(dá)、配體LFC等)來(lái)可視化輸出仁热，而不是使用circos圖可視化榜揖。無(wú)奈使用圈圖可視化的需求太多了，所以有了這篇教程抗蠢。

作為細(xì)胞相互作用的例子举哟，我們將使用Medaglia等人的小鼠NICHE-seq數(shù)據(jù)來(lái)探索淋巴細(xì)胞脈絡(luò)叢腦膜炎病毒(LCMV)感染前和72小時(shí)后腹股溝淋巴結(jié)T細(xì)胞區(qū)域的細(xì)胞間通訊[見(jiàn)@ medagli_spatial_2017]。我們將利用NicheNet來(lái)探索免疫細(xì)胞對(duì)這種LCMV感染的反應(yīng)迅矛。在本數(shù)據(jù)集中妨猩，觀察到穩(wěn)定狀態(tài)下的CD8 T細(xì)胞和LCMV感染后的CD8 T細(xì)胞之間的差異表達(dá)。NicheNet可以應(yīng)用于觀察淋巴結(jié)中的幾種免疫細(xì)胞群(如單核細(xì)胞秽褒、樹突狀細(xì)胞册赛、NK細(xì)胞、B細(xì)胞震嫉、CD4 T細(xì)胞)如何調(diào)節(jié)和誘導(dǎo)這些觀察到的基因表達(dá)變化森瘪。NicheNet將特別優(yōu)先考慮來(lái)自這些免疫細(xì)胞的配體及其在LCMV感染后表達(dá)改變的靶基因。

裝載所需的包票堵，用處理過(guò)的相互作用細(xì)胞的表達(dá)數(shù)據(jù)和NicheNet的配體-靶先驗(yàn)?zāi)Ｐ投蟛恰⑴潴w-受體網(wǎng)絡(luò)和加權(quán)綜合網(wǎng)絡(luò)讀取Seurat對(duì)象。

#devtools::install_github("saeyslab/nichenetr")
library(nichenetr)
library(Seurat)
library(tidyverse)
library(circlize)
ligand_target_matrix = readRDS(url("https://zenodo.org/record/3260758/files/ligand_target_matrix.rds"))
ligand_target_matrix[1:5,1:5]

                CXCL1        CXCL2        CXCL3        CXCL5         PPBP
A1BG     3.534343e-04 4.041324e-04 3.729920e-04 3.080640e-04 2.628388e-04
A1BG-AS1 1.650894e-04 1.509213e-04 1.583594e-04 1.317253e-04 1.231819e-04
A1CF     5.787175e-04 4.596295e-04 3.895907e-04 3.293275e-04 3.211944e-04
A2M      6.027058e-04 5.996617e-04 5.164365e-04 4.517236e-04 4.590521e-04
A2M-AS1  8.898724e-05 8.243341e-05 7.484018e-05 4.912514e-05 5.120439e-05

> lr_network = readRDS(url("https://zenodo.org/record/3260758/files/lr_network.rds"))
> head(lr_network)
# A tibble: 6 x 4
  from  to    source         database
  <chr> <chr> <chr>          <chr>   
1 CXCL1 CXCR2 kegg_cytokines kegg    
2 CXCL2 CXCR2 kegg_cytokines kegg    
3 CXCL3 CXCR2 kegg_cytokines kegg    
4 CXCL5 CXCR2 kegg_cytokines kegg    
5 PPBP  CXCR2 kegg_cytokines kegg    
6 CXCL6 CXCR2 kegg_cytokines kegg    

weighted_networks = readRDS(url("https://zenodo.org/record/3260758/files/weighted_networks.rds"))
head(weighted_networks$lr_sig) # interactions and their weights in the ligand-receptor + signaling network
# A tibble: 6 x 3
  from  to     weight
  <chr> <chr>   <dbl>
1 A1BG  ABCC6  0.422 
2 A1BG  ACE2   0.101 
3 A1BG  ADAM10 0.0970
4 A1BG  AGO1   0.0525
5 A1BG  AKT1   0.0855
6 A1BG  ANXA7  0.457 

seuratObj = readRDS(url("https://zenodo.org/record/3531889/files/seuratObj.rds"))
seuratObj@meta.data %>% head()
##         nGene nUMI orig.ident aggregate res.0.6 celltype nCount_RNA nFeature_RNA
## W380370   880 1611      LN_SS        SS       1    CD8 T       1607          876
## W380372   541  891      LN_SS        SS       0    CD4 T        885          536
## W380374   742 1229      LN_SS        SS       0    CD4 T       1223          737
## W380378   847 1546      LN_SS        SS       1    CD8 T       1537          838
## W380379   839 1606      LN_SS        SS       0    CD4 T       1603          836
## W380381   517  844      LN_SS        SS       0    CD4 T        840          513

seuratObj
An object of class Seurat 
13541 features across 5027 samples within 1 assay 
Active assay: RNA (13541 features, 1575 variable features)
 4 dimensional reductions calculated: cca, cca.aligned, tsne, pca

觀察存在哪些細(xì)胞群:CD4 T細(xì)胞(包括調(diào)節(jié)性T細(xì)胞)悴势、CD8 T細(xì)胞窗宇、B細(xì)胞、NK細(xì)胞特纤、樹突狀細(xì)胞(DCs)和炎癥單核細(xì)胞

seuratObj@meta.data$celltype %>% table() # note that the number of cells of some cell types is very low and should preferably be higher for a real application
## .
##     B CD4 T CD8 T    DC  Mono    NK  Treg 
##   382  2562  1645    18    90   131   199
DimPlot(seuratObj, reduction = "tsne",label = T)+ theme_bw()

seuratObj@meta.data$aggregate %>% table()
## .
## LCMV   SS 
## 3886 1141
DimPlot(seuratObj, reduction = "tsne", group.by = "aggregate")

對(duì)Seurat對(duì)象的NicheNet分析:解釋兩種條件的差異表達(dá)军俊。在這個(gè)案例研究中，接收細(xì)胞群是' CD8 T '細(xì)胞群捧存，而發(fā)送細(xì)胞群是' CD4 T '粪躬， ' Treg '， ' Mono '昔穴， ' NK '镰官， ' B '和' DC '。以上描述的功能將考慮一個(gè)基因在至少一個(gè)集群中預(yù)定義的部分細(xì)胞(默認(rèn)為10%)中表達(dá)時(shí)的表達(dá)吗货。我們感興趣的基因是LCMV感染后CD8 T細(xì)胞中差異表達(dá)的基因泳唠。因此，感興趣的條件是“LCMV”宙搬，而參考/穩(wěn)態(tài)條件是“SS”笨腥。條件的信息可以從元數(shù)據(jù)列“aggregate”中提取出來(lái)拓哺，計(jì)算差異基因的方法是標(biāo)準(zhǔn)的Seurat Wilcoxon檢驗(yàn)。

用于預(yù)測(cè)活性靶基因和構(gòu)建活性配體-受體網(wǎng)絡(luò)的配體的數(shù)量默認(rèn)為20個(gè)脖母。

# indicated cell types should be cell class identities
# check via: 
# seuratObj %>% Idents() %>% table()
sender_celltypes = c("CD4 T","Treg", "Mono", "NK", "B", "DC")
nichenet_output = nichenet_seuratobj_aggregate(
  seurat_obj = seuratObj, 
  receiver = "CD8 T", 
  condition_colname = "aggregate", condition_oi = "LCMV", condition_reference = "SS", 
  sender = sender_celltypes, 
  ligand_target_matrix = ligand_target_matrix, lr_network = lr_network, weighted_networks = weighted_networks, organism = "mouse")
## [1] "Read in and process NicheNet's networks"
## [1] "Define expressed ligands and receptors in receiver and sender cells"
## [1] "Perform DE analysis in receiver cell"
## [1] "Perform NicheNet ligand activity analysis"
## [1] "Infer active target genes of the prioritized ligands"
## [1] "Infer receptors of the prioritized ligands"

# 輸出的是一個(gè)列表：
 nichenet_output %>% names()
 [1] "ligand_activities"                "top_ligands"                      "top_targets"                     
 [4] "top_receptors"                    "ligand_target_matrix"             "ligand_target_heatmap"           
 [7] "ligand_target_df"                 "ligand_activity_target_heatmap"   "ligand_receptor_matrix"          
[10] "ligand_receptor_heatmap"          "ligand_receptor_df"               "ligand_receptor_matrix_bonafide" 
[13] "ligand_receptor_heatmap_bonafide" "ligand_receptor_df_bonafide"      "geneset_oi"                      
[16] "background_expressed_genes"      

配體活性分析結(jié)果士鸥。NicheNet做的第一件事，是根據(jù)預(yù)測(cè)的配體活性來(lái)確定配體的優(yōu)先級(jí)镶奉。使用如下命令查看這些配體的排名:

nichenet_output$ligand_activities
## # A tibble: 44 x 6
##    test_ligand auroc  aupr pearson  rank bona_fide_ligand
##    <chr>       <dbl> <dbl>   <dbl> <dbl> <lgl>           
##  1 Ebi3        0.662 0.238  0.219      1 FALSE           
##  2 Il15        0.596 0.160  0.109      2 TRUE            
##  3 Crlf2       0.560 0.160  0.0890     3 FALSE           
##  4 App         0.499 0.134  0.0750     4 TRUE            
##  5 Tgfb1       0.498 0.134  0.0631     5 TRUE            
##  6 Ptprc       0.539 0.142  0.0602     6 TRUE            
##  7 H2-M3       0.526 0.149  0.0533     7 TRUE            
##  8 Icam1       0.544 0.134  0.0496     8 TRUE            
##  9 Cxcl10      0.536 0.134  0.0457     9 TRUE            
## 10 Adam17      0.517 0.129  0.0378    10 TRUE            
## # ... with 34 more rows

這些配體由一個(gè)或多個(gè)輸入發(fā)送細(xì)胞表達(dá)。要看哪個(gè)細(xì)胞群表達(dá)了這些配體崭放，你可以運(yùn)行以下程序:

DotPlot(seuratObj, features = nichenet_output$top_ligands %>% rev(), cols = "RdYlBu") + RotatedAxis()

如你所見(jiàn)哨苛，大多數(shù)排名靠前的op配體似乎主要由樹突狀細(xì)胞和單核細(xì)胞表達(dá)。
觀察這些配體在LCMV感染后是否有差異表達(dá)也是很有趣的币砂。

DotPlot(seuratObj, features = nichenet_output$top_ligands %>% rev(), split.by = "aggregate") + RotatedAxis()

VlnPlot(seuratObj, features = c("Il15", "Cxcl10","Cxcl16"), split.by = "aggregate", pt.size = 0, combine = T)
## [[1]]

image.png

推斷活躍的配體-靶標(biāo)連接
NicheNet還推斷出這些頂級(jí)配體的活性靶基因建峭。要查看哪個(gè)頂級(jí)配體被預(yù)測(cè)調(diào)控了哪些差異表達(dá)基因的表達(dá)，可以運(yùn)行以下命令來(lái)查看熱圖:

nichenet_output$ligand_target_heatmap

Circos繪圖來(lái)可視化配體-靶標(biāo)和配體-受體的相互作用决摧。這一可視化分組根據(jù)最強(qiáng)表達(dá)的細(xì)胞類型預(yù)測(cè)活性配體亿蒸。因此，我們需要確定每種細(xì)胞類型掌桩，它們表達(dá)的配體比其他細(xì)胞類型更強(qiáng)边锁。計(jì)算發(fā)送細(xì)胞中平均配體表達(dá)量。

# avg_expression_ligands = AverageExpression(seuratObj %>% subset(subset = aggregate == "LCMV"),features = nichenet_output$top_ligands) # if want to look specifically in LCMV-only cells
avg_expression_ligands = AverageExpression(seuratObj, features = nichenet_output$top_ligands)

分配配體給發(fā)送細(xì)胞
為了給發(fā)送端細(xì)胞類型分配配體波岛，我們可以查找哪個(gè)發(fā)送端細(xì)胞類型的表達(dá)式高于平均值+ SD茅坛。

sender_ligand_assignment = avg_expression_ligands$RNA %>% apply(1, function(ligand_expression){
  ligand_expression > (ligand_expression %>% mean() + ligand_expression %>% sd())
  }) %>% t()

sender_ligand_assignment[1:4,1:4]
      CD8 T CD4 T  Treg     B
Ebi3  FALSE FALSE FALSE FALSE
Il15  FALSE FALSE FALSE FALSE
Crlf2 FALSE FALSE FALSE FALSE
App   FALSE FALSE FALSE FALSE

sender_ligand_assignment = sender_ligand_assignment %>% apply(2, function(x){x[x == TRUE]}) %>% purrr::keep(function(x){length(x) > 0})
names(sender_ligand_assignment)
## [1] "B"    "NK"   "Mono" "DC"

(sender_ligand_assignment)
$B
H2-M3 
 TRUE 

$NK
Ptprc Itgb1 
 TRUE  TRUE 

$Mono
  Ebi3  Crlf2    App  Tgfb1 Cxcl10 Adam17 Cxcl11  Cxcl9 Sema4d     C3  Anxa1 
  TRUE   TRUE   TRUE   TRUE   TRUE   TRUE   TRUE   TRUE   TRUE   TRUE   TRUE 

$DC
  Il15  Icam1 H2-T23   Ccl5 Cxcl16  Itgb1 
  TRUE   TRUE   TRUE   TRUE   TRUE   TRUE 

頂部的配體似乎在B細(xì)胞、NK細(xì)胞则拷、單核細(xì)胞和DCs中表達(dá)最強(qiáng)烈贡蓖。我們也會(huì)知道在多種細(xì)胞類型中哪些配體是常見(jiàn)的(=特定于> 1細(xì)胞類型的配體，或在之前的代碼塊中未指定給某個(gè)細(xì)胞類型的配體)』筒纾現(xiàn)在確定哪些優(yōu)先配體是由CAFs或內(nèi)皮細(xì)胞表達(dá)的斥铺。

all_assigned_ligands = sender_ligand_assignment %>% lapply(function(x){names(x)}) %>% unlist()
unique_ligands = all_assigned_ligands %>% table() %>% .[. == 1] %>% names()
general_ligands = nichenet_output$top_ligands %>% setdiff(unique_ligands)

B_specific_ligands = sender_ligand_assignment$B %>% names() %>% setdiff(general_ligands)
NK_specific_ligands = sender_ligand_assignment$NK %>% names() %>% setdiff(general_ligands)
Mono_specific_ligands = sender_ligand_assignment$Mono %>% names() %>% setdiff(general_ligands)
DC_specific_ligands = sender_ligand_assignment$DC %>% names() %>% setdiff(general_ligands)

ligand_type_indication_df = tibble(
  ligand_type = c(rep("B-specific", times = B_specific_ligands %>% length()),
                  rep("NK-specific", times = NK_specific_ligands %>% length()),
                  rep("Mono-specific", times = Mono_specific_ligands %>% length()),
                  rep("DC-specific", times = DC_specific_ligands %>% length()),
                  rep("General", times = general_ligands %>% length())),
  ligand = c(B_specific_ligands, NK_specific_ligands, Mono_specific_ligands, DC_specific_ligands, general_ligands))

ligand_type_indication_df %>% head
# A tibble: 6 x 2
  ligand_type   ligand
  <chr>         <chr> 
1 B-specific    H2-M3 
2 NK-specific   Ptprc 
3 Mono-specific Ebi3  
4 Mono-specific Crlf2 
5 Mono-specific App   
6 Mono-specific Tgfb1 

定義感興趣的配體-目標(biāo)鏈接
為了避免circos圖中有太多配體目標(biāo)鏈接，我們將只顯示權(quán)重高于預(yù)定義截止值的鏈接:屬于最低分?jǐn)?shù)的40%的鏈接被刪除坛善。這并不是說(shuō)用于這種可視化的邊界和其他邊界可以根據(jù)用戶的需要進(jìn)行更改晾蜘。

active_ligand_target_links_df = nichenet_output$ligand_target_df %>% mutate(target_type = "LCMV-DE") %>% inner_join(ligand_type_indication_df) # if you want ot make circos plots for multiple gene sets, combine the different data frames and differentiate which target belongs to which gene set via the target type

cutoff_include_all_ligands = active_ligand_target_links_df$weight %>% quantile(0.40)

active_ligand_target_links_df_circos = active_ligand_target_links_df %>% filter(weight > cutoff_include_all_ligands)

ligands_to_remove = setdiff(active_ligand_target_links_df$ligand %>% unique(), active_ligand_target_links_df_circos$ligand %>% unique())
targets_to_remove = setdiff(active_ligand_target_links_df$target %>% unique(), active_ligand_target_links_df_circos$target %>% unique())
  
circos_links = active_ligand_target_links_df %>% filter(!target %in% targets_to_remove &!ligand %in% ligands_to_remove)

 circos_links
# A tibble: 124 x 5
   ligand target  weight target_type ligand_type  
   <chr>  <chr>    <dbl> <chr>       <chr>        
 1 Ebi3   Cd274  0.00325 LCMV-DE     Mono-specific
 2 Ebi3   Cd53   0.00321 LCMV-DE     Mono-specific
 3 Ebi3   Ddit4  0.00335 LCMV-DE     Mono-specific
 4 Ebi3   Id3    0.00373 LCMV-DE     Mono-specific
 5 Ebi3   Ifit3  0.00320 LCMV-DE     Mono-specific
 6 Ebi3   Irf1   0.00692 LCMV-DE     Mono-specific
 7 Ebi3   Irf7   0.00312 LCMV-DE     Mono-specific
 8 Ebi3   Irf9   0.00543 LCMV-DE     Mono-specific
 9 Ebi3   Parp14 0.00336 LCMV-DE     Mono-specific
10 Ebi3   Pdcd4  0.00335 LCMV-DE     Mono-specific
# ... with 114 more rows

準(zhǔn)備circos可視化:給每個(gè)片段配體和目標(biāo)特定的顏色和順序

grid_col_ligand =c("General" = "lawngreen",
            "NK-specific" = "royalblue",
            "B-specific" = "darkgreen",
            "Mono-specific" = "violet",
            "DC-specific" = "steelblue2")
grid_col_target =c(
            "LCMV-DE" = "tomato")

grid_col_tbl_ligand = tibble(ligand_type = grid_col_ligand %>% names(), color_ligand_type = grid_col_ligand)
grid_col_tbl_target = tibble(target_type = grid_col_target %>% names(), color_target_type = grid_col_target)

circos_links = circos_links %>% mutate(ligand = paste(ligand," ")) # extra space: make a difference between a gene as ligand and a gene as target!
circos_links = circos_links %>% inner_join(grid_col_tbl_ligand) %>% inner_join(grid_col_tbl_target)
links_circle = circos_links %>% select(ligand,target, weight)

ligand_color = circos_links %>% distinct(ligand,color_ligand_type)
grid_ligand_color = ligand_color$color_ligand_type %>% set_names(ligand_color$ligand)
target_color = circos_links %>% distinct(target,color_target_type)
grid_target_color = target_color$color_target_type %>% set_names(target_color$target)

grid_col =c(grid_ligand_color,grid_target_color)

# give the option that links in the circos plot will be transparant ~ ligand-target potential score
transparency = circos_links %>% mutate(weight =(weight-min(weight))/(max(weight)-min(weight))) %>% mutate(transparency = 1-weight) %>% .$transparency 

準(zhǔn)備可視化的circos:排序配體和目標(biāo)

target_order = circos_links$target %>% unique()
ligand_order = c(Mono_specific_ligands, DC_specific_ligands, NK_specific_ligands,B_specific_ligands, general_ligands) %>% c(paste(.," ")) %>% intersect(circos_links$ligand)
order = c(ligand_order,target_order)

準(zhǔn)備circos可視化:定義不同片段之間的間隙

width_same_cell_same_ligand_type = 0.5
width_different_cell = 6
width_ligand_target = 15
width_same_cell_same_target_type = 0.5

gaps = c(
  # width_ligand_target,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "Mono-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "DC-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "NK-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "B-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow() -1)),
  width_ligand_target,
  rep(width_same_cell_same_target_type, times = (circos_links %>% filter(target_type == "LCMV-DE") %>% distinct(target) %>% nrow() -1)),
  width_ligand_target
  )

渲染circos的情節(jié)(所有鏈接相同的透明度)。只有表明每個(gè)靶基因的阻滯的寬度與配體-靶的調(diào)控電位成正比(~支持調(diào)控相互作用的先驗(yàn)知識(shí))眠屎。

circos.par(gap.degree = gaps)
chordDiagram(links_circle, directional = 1,order=order,link.sort = TRUE, link.decreasing = FALSE, grid.col = grid_col,transparency = 0, diffHeight = 0.005, direction.type = c("diffHeight", "arrows"),link.arr.type = "big.arrow", link.visible = links_circle$weight >= cutoff_include_all_ligands,annotationTrack = "grid", 
    preAllocateTracks = list(track.height = 0.075))
# we go back to the first track and customize sector labels
circos.track(track.index = 1, panel.fun = function(x, y) {
    circos.text(CELL_META$xcenter, CELL_META$ylim[1], CELL_META$sector.index,
        facing = "clockwise", niceFacing = TRUE, adj = c(0, 0.55), cex = 1)
}, bg.border = NA) #

circos.clear()

繪制circos圖(透明度由配體-靶標(biāo)相互作用的調(diào)控潛力決定)

circos.par(gap.degree = gaps)
chordDiagram(links_circle, directional = 1,order=order,link.sort = TRUE, link.decreasing = FALSE, grid.col = grid_col,transparency = transparency, diffHeight = 0.005, direction.type = c("diffHeight", "arrows"),link.arr.type = "big.arrow", link.visible = links_circle$weight >= cutoff_include_all_ligands,annotationTrack = "grid", 
    preAllocateTracks = list(track.height = 0.075))
# we go back to the first track and customize sector labels
circos.track(track.index = 1, panel.fun = function(x, y) {
    circos.text(CELL_META$xcenter, CELL_META$ylim[1], CELL_META$sector.index,
        facing = "clockwise", niceFacing = TRUE, adj = c(0, 0.55), cex = 1)
}, bg.border = NA) #

circos.clear()

svg("ligand_target_circos.svg", width = 10, height = 10)
circos.par(gap.degree = gaps)
chordDiagram(links_circle, directional = 1,order=order,link.sort = TRUE, link.decreasing = FALSE, grid.col = grid_col,transparency = transparency, diffHeight = 0.005, direction.type = c("diffHeight", "arrows"),link.arr.type = "big.arrow", link.visible = links_circle$weight >= cutoff_include_all_ligands,annotationTrack = "grid",
    preAllocateTracks = list(track.height = 0.075))
# we go back to the first track and customize sector labels
circos.track(track.index = 1, panel.fun = function(x, y) {
    circos.text(CELL_META$xcenter, CELL_META$ylim[1], CELL_META$sector.index,
        facing = "clockwise", niceFacing = TRUE, adj = c(0, 0.55), cex = 1)
}, bg.border = NA) #
circos.clear()
dev.off()

在circos圖中可視化優(yōu)先配體與受體的相互作用

lr_network_top_df = nichenet_output$ligand_receptor_df %>% mutate(receptor_type = "LCMV_CD8T_receptor") %>% inner_join(ligand_type_indication_df)
grid_col_ligand =c("General" = "lawngreen",
            "NK-specific" = "royalblue",
            "B-specific" = "darkgreen",
            "Mono-specific" = "violet",
            "DC-specific" = "steelblue2")
grid_col_receptor =c(
            "LCMV_CD8T_receptor" = "darkred")

grid_col_tbl_ligand = tibble(ligand_type = grid_col_ligand %>% names(), color_ligand_type = grid_col_ligand)
grid_col_tbl_receptor = tibble(receptor_type = grid_col_receptor %>% names(), color_receptor_type = grid_col_receptor)

circos_links = lr_network_top_df %>% mutate(ligand = paste(ligand," ")) # extra space: make a difference between a gene as ligand and a gene as receptor!
circos_links = circos_links %>% inner_join(grid_col_tbl_ligand) %>% inner_join(grid_col_tbl_receptor)
links_circle = circos_links %>% select(ligand,receptor, weight)

ligand_color = circos_links %>% distinct(ligand,color_ligand_type)
grid_ligand_color = ligand_color$color_ligand_type %>% set_names(ligand_color$ligand)
receptor_color = circos_links %>% distinct(receptor,color_receptor_type)
grid_receptor_color = receptor_color$color_receptor_type %>% set_names(receptor_color$receptor)

grid_col =c(grid_ligand_color,grid_receptor_color)

# give the option that links in the circos plot will be transparant ~ ligand-receptor potential score
transparency = circos_links %>% mutate(weight =(weight-min(weight))/(max(weight)-min(weight))) %>% mutate(transparency = 1-weight) %>% .$transparency 

制備可視化的circos:有序配體和受體

receptor_order = circos_links$receptor %>% unique()
ligand_order = c(Mono_specific_ligands, DC_specific_ligands, NK_specific_ligands,B_specific_ligands, general_ligands) %>% c(paste(.," ")) %>% intersect(circos_links$ligand)
order = c(ligand_order,receptor_order)

準(zhǔn)備馬戲團(tuán)可視化:定義不同片段之間的間隙

width_same_cell_same_ligand_type = 0.5
width_different_cell = 6
width_ligand_receptor = 15
width_same_cell_same_receptor_type = 0.5

gaps = c(
  # width_ligand_target,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "Mono-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "DC-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "NK-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "B-specific") %>% distinct(ligand) %>% nrow() -1)),
  width_different_cell,
  rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow() -1)),
  width_ligand_receptor,
  rep(width_same_cell_same_receptor_type, times = (circos_links %>% filter(receptor_type == "LCMV_CD8T_receptor") %>% distinct(receptor) %>% nrow() -1)),
  width_ligand_receptor
  )

渲染馬戲團(tuán)的情節(jié)(所有鏈接相同的透明度)笙纤。只有表明每個(gè)受體的阻滯的寬度與配體-受體相互作用的重量成比例(~支持相互作用的先驗(yàn)知識(shí))。

circos.par(gap.degree = gaps)
chordDiagram(links_circle, directional = 1,order=order,link.sort = TRUE, link.decreasing = FALSE, grid.col = grid_col,transparency = 0, diffHeight = 0.005, direction.type = c("diffHeight", "arrows"),link.arr.type = "big.arrow", link.visible = links_circle$weight >= cutoff_include_all_ligands,annotationTrack = "grid", 
    preAllocateTracks = list(track.height = 0.075))
# we go back to the first track and customize sector labels
circos.track(track.index = 1, panel.fun = function(x, y) {
    circos.text(CELL_META$xcenter, CELL_META$ylim[1], CELL_META$sector.index,
        facing = "clockwise", niceFacing = TRUE, adj = c(0, 0.55), cex = 0.8)
}, bg.border = NA) #

circos.clear()

渲染circos圖(透明程度由配體-受體相互作用的先驗(yàn)相互作用權(quán)重決定——正如指示每個(gè)受體的塊的寬度)

circos.par(gap.degree = gaps)
chordDiagram(links_circle, directional = 1,order=order,link.sort = TRUE, link.decreasing = FALSE, grid.col = grid_col,transparency = transparency, diffHeight = 0.005, direction.type = c("diffHeight", "arrows"),link.arr.type = "big.arrow", link.visible = links_circle$weight >= cutoff_include_all_ligands,annotationTrack = "grid", 
    preAllocateTracks = list(track.height = 0.075))
# we go back to the first track and customize sector labels
circos.track(track.index = 1, panel.fun = function(x, y) {
    circos.text(CELL_META$xcenter, CELL_META$ylim[1], CELL_META$sector.index,
        facing = "clockwise", niceFacing = TRUE, adj = c(0, 0.55), cex = 0.8)
}, bg.border = NA) #

circos.clear()

Perform NicheNet analysis starting from a Seurat object:vignette("seurat_wrapper", package="nichenetr")
Circos plot visualization to show active ligand-target links between interacting cells:
https://github.com/saeyslab/nichenetr/blob/master/vignettes/seurat_wrapper_circos.md
https://github.com/saeyslab/nichenetr/issues/5