Circos plot for Seurat single cell tutorial?

[Seandelao]

Hello,

Thanks again for the awesome package, I am finding some really cool things with it!

I was wondering if you were planning on adding a circos plot tutorial for the Seurat single cell vignette? I was trying to follow the current example utilizing my Seurat NicheNetR analysis, but I think I am a little lost on which of the NicheNet outputs to use for the different steps.

I think an example vignette would be super helpful! Thanks again!

[browaeysrobin]

Hi @Seandelao,

We will be working on this in the near future.

[saeedfc]

@browaeysrobin ,

Thanks for he excellent tool. Another suggestion on the same line will be to include in the tutorial, an option to have additional section layer in the circos plot to indicate the receptor corresponding to the target genes as you showed in the 6B figure of the article 'Stellate Cells, Hepatocytes, and Endothelial Cells Imprint the Kupffer Cell Identity on Monocytes Colonizing the Liver Macrophage Niche' from your lab. It might violate copyrights if I past the image here. Sorry for that.

[ccruizm]

Hello @browaeysrobin!

I have been playing around with your nice tool! The tutorials have been a great help. I am also trying to generate the circos plot based on a Seurat object and 10x data. I am stuck now on how to prioritize ligands expressed by specific cell types. In your vignette, you first filtered them based on their log2(average(TPM(i)1…k)+1) value (>4) and later, by a difference of 2 expression units (for prioritization). How would you suggest to do it with data from 10x data stored in a Seurat object? In the beginning, one filter considering a gene to be expressed in at least 10% of cells in one cluster, but in this particular step, It is not clear to me how to prioritize the genes for downstream analysis.

Thanks in advance for your help.

[browaeysrobin]

Hi @ccruizm

These steps are not really for prioritization, but just to assign ligands to one of the possible sender cell types (i.e. the sender cell type that most strongly expresses this ligand) or to the group of 'generally' expressed ligands.
You could calculate the average expression of the ligand in each cell type / cluster. You can write a function yourself for this but there is also the option to use the AverageExpression function from the Seurat package. Then you can find the cell-type specific ligands by looking at whether the average expression of a ligand in cluster is X higher than in other clusters. For the shown data, we used a difference of 2 units, but here you should probably use a less stringent cutoff. (Or alternatively, you could just assign a ligand to its highest expressing cluster)

[ccruizm]

Thanks for the suggestions @browaeysrobin! I followed your advice and seems to work. Great package!

[rsggsr]

Hi, @browaeysrobin

Thanks for this guideline! Took me a while for figuring out how to do it but it finally works! Thanks for this interesting and powerful tool!!

[paolo-kunderfranco]

Hi all, fantastic package!
Any update on the circos plot tutorial for the Seurat single cell vignette?
Best
Paolo

[yong114]

Thanks for the suggestions @browaeysrobin! I followed your advice and seems to work. Great package!

Hi @ccruizm @rsggsr

Would you mind to share the codes please for generating the Circos plot based on Seurat object? It would be very helpful if you could share those please.

Thank you in advance and much appreciated. (Thanks @browaeysrobin for this awesome tool!)

[Elo-mars]

Thanks for the suggestions @browaeysrobin! I followed your advice and seems to work. Great package!

Hello @ccruizm,
like others, I'm not sure I understand what @browaeysrobin explained in his post ... (sorry)
Do you have to specify one ligand of preference? I would like to remained unbiased but then not sure how I can apply AverageExpression ...
Would you mind sharing your code to design circle plots from a Seurat Object?

I have this, so I have the most active ligands, their sender cell pop, the receptors and the target genes modified by all this, it would be nice to also have a circle plot to summarize all that

Thanks!
Elo
Combined_plots_3.pdf

[browaeysrobin]

Hello @Elo-mars, @paolo-kunderfranco , @TrumanZYX, @yong114

I quickly made an example vignette that combines the Seurat wrapper with the circos plot visualization.
https://github.com/saeyslab/nichenetr/blob/master/vignettes/seurat_wrapper_circos.md

Note that I generally recommend a visualization of the output by combining several heatmaps (ligand activity, ligand-target links, ligand-receptor links, ligand expression, ligand LFC,…) over using a circos plot visualization. Certainly for cases with many sender cell types and ligands that are expressed by more than one sender cell type. Because in those cases, the circos plot is much less informative and could lead to wrong interpretation of the results.

[Elo-mars]

Super thanks Robin! I'll try that this week :)

[rencav18]

Hey @browaeysrobin

First and foremost, I want to thank you for the Circos Seurat Plot vignette! I'm relatively new to coding and everything so I appreciate it very much, I know it takes a lot of hard work, time, and effort, so thank you again!

For the Circos plot Seurat, I seem to be getting an error in my script, most specifically with the gaps argument. For the section with the header "Prepare the circos visualization: define the gaps between the different segments" I keep getting
"Error in rep(width_same_cell_same_ligand_type, times = (circos_links %>% :
invalid 'times' argument" Do you know of any way to fix this? Do you need to see the script? (Sorry new to coding)

Separate question, using this Circos Plot Seurat Vignette, is it ever possible to compare between Single cell RNA sequencing data set and bulk RNA sequencing data set by normalizing the average expression values from one to the other?

Once again,
Thank you!

[browaeysrobin]

Hi @rencav18

Your error: probably a coding mistake by you. The error message says that you have an invalid 'times' argument, so check whether what you have put after the times argument is working code and returns a number.

Could you clarify your question about comparing scRNAseq and bulk data? I don't see what you try to get out of this in the NicheNet context.

[whitneyt1]

@rencav18

I had the same error and ran unique(circos_links$ligand_type) and limited my gaps commands to only include those and it worked!

[DRSEI]

hello guys
I need to reopen this case again: I am trying to do circos plot and I am getting a 2 type of error:

I have attached my first error in below :

active_ligand_target_links_df = nichenet_output$ligand_target_df %>% mutate(target_type = "Day29PS-DE") %>% inner_join(ligand_type_indication_df)

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

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)

##Prepare the circos visualization: give each segment of ligands and targets a specific color and order
grid_col_ligand =c("General" = "lawngreen",

•                "pDC-specific" = "royalblue",
  

•                "cDC-specific" = "darkgreen",
  

•                "Mono-specific" = "violet",
  

•                "CD16-specific" = "steelblue2")
  

grid_col_target =c(

• "Day29PS-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)
Joining, by = "ligand_type"
Joining, by = "target_type"
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

target_order = circos_links$target %>% unique()
ligand_order = c(pDC_specific_ligands, cDC_specific_ligands, Mono_specific_ligands, CD16_specific_ligands, general_ligands) %>% c(paste(.," ")) %>% intersect(circos_links$ligand)
order = c(ligand_order,target_order)

width_same_cell_same_ligand_type = 0.5
width_different_cell = 6
width_ligand_target = 15
width_same_cell_same_target_type = 0.5
unique(circos_links$ligand_type)
[1] "General" "Mono-specific" "CD16-specific" "cDC-specific" "pDC-specific"
gaps = c(

• width_ligand_target,

• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• 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 == "CD16-specific") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "cDC-specific") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "pDC-specific" ) %>% distinct(ligand) %>% nrow() -1)),
• width_ligand_target,
• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "Day29PS-DE" ) %>% distinct(ligand) %>% nrow() -1)),
• width_ligand_target
• )
  **Error in rep(width_same_cell_same_ligand_type, times = (circos_links %>% :
  invalid 'times' argument

**

When I removed the ( width_ligand_target,

• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "Day29PS-DE" ) %>% distinct(ligand) %>% nrow() -1)),

then i am getting error bellow :

unique(circos_links$ligand_type)
[1] "General" "Mono-specific" "CD16-specific" "cDC-specific" "pDC-specific"
gaps = c(

• width_ligand_target,

• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• 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 == "CD16-specific") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "cDC-specific") %>% distinct(ligand) %>% nrow() -1)),
• width_different_cell,
• rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "pDC-specific" ) %>% distinct(ligand) %>% nrow() -1)),
• width_ligand_target
• )

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))
  

Error: Since gap.degree parameter has length larger than 1, it should have same length as the number of sectors.

Thank you in advance and much appreciated

[browaeysrobin]

Dear @DRSEI

The problem is in 1 or more of the times argument of the rep function to define the gaps variable.
I suggest you print out each of these times arguments.

eg circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow()
This should return a number. If not check the output of circos_links %>% filter(ligand_type == "General") to see whether this is not an empty dataframe

Do this for every ligand_type you have

[DRSEI]

Dear @browaeysrobin

I have followed your instruction and Still getting below error
Since gap.degree parameter has length larger than 1, it should have same length as the number of sectors.

> circos_links %>% filter(ligand_type == "General") %>% distinct(ligand) %>% nrow()
[1] 11
> circos_links %>% filter(ligand_type == "Mono-specific") %>% distinct(ligand) %>% nrow()
[1] 1
> circos_links %>% filter(ligand_type == "CD16-specific") %>% distinct(ligand) %>% nrow()
[1] 2
> circos_links %>% filter(ligand_type == "cDC-specific") %>% distinct(ligand) %>% nrow()
[1] 1
> circos_links %>% filter(ligand_type == "pDC-specific") %>% distinct(ligand) %>% nrow()
[1] 1


> circos_links %>% filter(ligand_type == "General")
# A tibble: 74 x 7
   ligand    target   weight target_type ligand_type color_ligand_type color_target_type
   <chr>     <chr>     <dbl> <chr>       <chr>       <chr>             <chr>            
 1 "TGFB1  " CCL3    0.00676 Day28-DE    General     lawngreen         tomato           
 2 "TGFB1  " GZMB    0.00595 Day28-DE    General     lawngreen         tomato           
 3 "TGFB1  " IFIT3   0.00467 Day28-DE    General     lawngreen         tomato           
 4 "TGFB1  " ITGB1   0.00466 Day28-DE    General     lawngreen         tomato           
 5 "TGFB1  " JUNB    0.00574 Day28-DE    General     lawngreen         tomato           
 6 "TGFB1  " KLF6    0.00481 Day28-DE    General     lawngreen         tomato           
 7 "IL1B  "  CCL3    0.00530 Day28-DE    General     lawngreen         tomato           
 8 "IL1B  "  GADD45B 0.00512 Day28-DE    General     lawngreen         tomato           
 9 "IL1B  "  ITGB1   0.00596 Day28-DE    General     lawngreen         tomato           
10 "IL1B  "  NFKBIA  0.00480 Day28-DE    General     lawngreen         tomato           
# ... with 64 more rows
> circos_links %>% filter(ligand_type == "Mono-specific")
# A tibble: 8 x 7
  ligand    target    weight target_type ligand_type   color_ligand_type color_target_type
  <chr>     <chr>      <dbl> <chr>       <chr>         <chr>             <chr>            
1 "CXCL2  " ATF3    0.00111  Day28-DE    Mono-specific violet            tomato           
2 "CXCL2  " DUSP1   0.00115  Day28-DE    Mono-specific violet            tomato           
3 "CXCL2  " GADD45B 0.00127  Day28-DE    Mono-specific violet            tomato           
4 "CXCL2  " JUNB    0.00128  Day28-DE    Mono-specific violet            tomato           
5 "CXCL2  " KLF6    0.000986 Day28-DE    Mono-specific violet            tomato           
6 "CXCL2  " NFKBIA  0.00132  Day28-DE    Mono-specific violet            tomato           
7 "CXCL2  " PFN1    0.00102  Day28-DE    Mono-specific violet            tomato           
8 "CXCL2  " TNFAIP3 0.00108  Day28-DE    Mono-specific violet            tomato           
> circos_links %>% filter(ligand_type == "CD16-specific")
# A tibble: 8 x 7
  ligand    target   weight target_type ligand_type   color_ligand_type color_target_type
  <chr>     <chr>     <dbl> <chr>       <chr>         <chr>             <chr>            
1 "IL15  "  GADD45B 0.00443 Day28-DE    CD16-specific steelblue2        tomato           
2 "IL15  "  GZMB    0.00453 Day28-DE    CD16-specific steelblue2        tomato           
3 "IL15  "  JUNB    0.00397 Day28-DE    CD16-specific steelblue2        tomato           
4 "CLCF1  " ATF3    0.00177 Day28-DE    CD16-specific steelblue2        tomato           
5 "CLCF1  " DUSP1   0.00160 Day28-DE    CD16-specific steelblue2        tomato           
6 "CLCF1  " GADD45B 0.00201 Day28-DE    CD16-specific steelblue2        tomato           
7 "CLCF1  " JUNB    0.00212 Day28-DE    CD16-specific steelblue2        tomato           
8 "CLCF1  " KLF6    0.00161 Day28-DE    CD16-specific steelblue2        tomato           
> circos_links %>% filter(ligand_type == "cDC-specific")
# A tibble: 7 x 7
  ligand    target   weight target_type ligand_type  color_ligand_type color_target_type
  <chr>     <chr>     <dbl> <chr>       <chr>        <chr>             <chr>            
1 "ALCAM  " DUSP1   0.00177 Day28-DE    cDC-specific darkgreen         tomato           
2 "ALCAM  " GADD45B 0.00217 Day28-DE    cDC-specific darkgreen         tomato           
3 "ALCAM  " JUNB    0.00206 Day28-DE    cDC-specific darkgreen         tomato           
4 "ALCAM  " KLF6    0.00158 Day28-DE    cDC-specific darkgreen         tomato           
5 "ALCAM  " NFKBIA  0.00215 Day28-DE    cDC-specific darkgreen         tomato           
6 "ALCAM  " PFN1    0.00203 Day28-DE    cDC-specific darkgreen         tomato           
7 "ALCAM  " TNFAIP3 0.00210 Day28-DE    cDC-specific darkgreen         tomato           
> circos_links %>% filter(ligand_type == "pDC-specific")
# A tibble: 8 x 7
  ligand   target    weight target_type ligand_type  color_ligand_type color_target_type
  <chr>    <chr>      <dbl> <chr>       <chr>        <chr>             <chr>            
1 "SELL  " ATF3    0.00124  Day28-DE    pDC-specific royalblue         tomato           
2 "SELL  " DUSP1   0.00110  Day28-DE    pDC-specific royalblue         tomato           
3 "SELL  " GADD45B 0.00118  Day28-DE    pDC-specific royalblue         tomato           
4 "SELL  " JUNB    0.00125  Day28-DE    pDC-specific royalblue         tomato           
5 "SELL  " KLF6    0.00101  Day28-DE    pDC-specific royalblue         tomato           
6 "SELL  " NFKBIA  0.00117  Day28-DE    pDC-specific royalblue         tomato           
7 "SELL  " PFN1    0.000944 Day28-DE    pDC-specific royalblue         tomato           
8 "SELL  " TNFAIP3 0.000944 Day28-DE    pDC-specific royalblue         tomato

> 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))

Error: Since gap.degree parameter has length larger than 1, it should have same length as the number of sectors.

[browaeysrobin]

Dear @DRSEI

I think this is the issue:
You seem to have only one ligand for some cell types. This will give some problems in defining the gaps as described in the vignette.
Cf:
rep(width_same_cell_same_ligand_type, times = (circos_links %>% filter(ligand_type == "Mono-specific") %>% distinct(ligand) %>% nrow() -1))
You see that the times argument is the nr of unique ligands subtracted by 1; which will give 0 if you have only one ligand for that celltype. 0 is an invalid argument of times in the rep function used. This makes sense of course, since this line of code is to make gaps between the sectors of the ligands in the same cell type. If you have only one ligand for a celltype, you should not define these gaps. Therefore I suggest you change the gaps vector in such a way that you remove lines of code that encode the gaps between ligands of one cell type in case you only have one ligand.

If that would not be the issue:

Can you check whether one of the ligands of one cell type is also a ligand in another cell type?
If yes, you should remove or add a space to that/those ligand(s).
I will be working in the near future to avoid this type of bug.

If that would still not be the issue:

Then I don't directly see what went wrong here. I would suggest that you try to locate the problem more specifically: check the number of sectors, the gaps vector etc.

[cenk-celik]

I encountered the same problem that @DRSEI reported here. I was able to draw the circos plot by simple not running circos.par(gap.degree = gaps). Basically, just generate circos plot without defining circos.par.

chordDiagram(links_circle, directional = 1, order = order, link.sort = T, link.decreasing = F, 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))
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 = T, adj = c(0, 0.55), cex = 0.5)
}, bg.border = NA)
circos.clear()