基于分类变量的集群圆形热图

Question

我想绘制一个类似于下图的圆形热图（从这里获得）。但是，我想按

clin.info$subtype

列进行聚类，该列与

meth

和

exp

的列名称相对应。

library(circlize)
library(ComplexHeatmap)

# Split by subtype category
split <- clin.info$subtype

# Colors
col.meth <-colorRamp2(c(-2, 0, 2), c("#00539C", "#f1f1f1", "#EEA47F"))
col.expr <- colorRamp2(c(-2, 0, 2), c("#317773", "#f1f1f1", "#E2D1F9"))
col.dnam.effect <- c("Attenuating"="#EC449B", "Enhancing"="#99F443", "ns"="#f1f1f1")
col.tf.role <- c("Repressor"="#CCF381", "Activator"="#4831D4", "ns"="#f1f1f1")
col.dnam.pval <- colorRamp2(c(-2, 0, 2), c("#f1f1f1", "#ffa3cf", "#FF69B4")) # for results.stratified.model
col.TF.pval <- colorRamp2(c(-2, 0, 2), c("#f1f1f1", "#EED971", "#FFA351")) # for results.stratified.model
col.synergistic.dnam.tf.pval <- colorRamp2(c(-2, 0, 2), c("#8A307F", "white", "#6883BC")) 
col.dist <- colorRamp2(c(0, 5000, 10000), c("#00FFFF", "#b1eeec", "#f1f1f1"))

# Plot function
circlize_plot = function() {
  chordDiagram(meth, split=split, col=col.meth, track.height=0.12)
  chordDiagram(distal.dmr.cpgs$DNAm.effect, col=col.dnam.effect, track.height=0.01)
  chordDiagram(distal.dmr.cpgs$TF.role, col=col.tf.role, track.height=0.01)
  chordDiagram(exp, col=col_expr, track.height=0.12)
  chordDiagram(distal.dmr.cpgs$RLM_DNAmGroup_pvalue, col=col.dnam.pval, track.height=0.01)
  chordDiagram(distal.dmr.cpgs$RLM_TF_pvalue, col=col.TF.pval, track.height=0.01)
  chordDiagram(distal.dmr.cpgs$`RLM_DNAmGroup:TF_pvalue`, col=col.synergistic.dnam.tf.pval, track.height=0.01)
  chordDiagram(distal.dmr.cpgs$distance_region_target_tss, col=col.dist, track.height = 0.01)
  circos.clear()
}

# Legends
lgd.subtype <- Legend(title="Subtype", legend_gp=gpar(fill=clust.col), at=names(clust.col))
lgd.meth <- Legend(title="Methylation", col_fun=col.meth)
lgd.expr <- Legend(title="Expression", col_fun=col.expr)
lgd.dnam.effect <- Legend(title="DNAm effect", at=names(col.dnam.effect), legend_gp=gpar(fill=col.dnam.effect))
lgd.tf.role <- Legend(title="TF role", at=names(col.tf.role), legend_gp=gpar(fill=col.tf.role))
lgd.dnam.pval <- Legend(title="DNAm P-value", col_fun=col.dnam.pval, at=c(0, 2, 4), labels=c(1, 0.01, 0.0001))
lgd.TF.pval <- Legend(title="TF P-value", col_fun=col.TF.pval, at=c(0, 2, 4), labels=c(1, 0.01, 0.0001))
lgd.synergistic.dnam.tf.pval <- Legend(title="Synergistic DNAm-TF P-value", col_fun=col.synergistic.dnam.tf.pval, at=c(0, 2, 4), labels=c(1, 0.01, 0.0001))
lgd.dist <- Legend(title="Dist to TSS", col_fun=col.dist, at=c(0, 5000, 10000), labels=c("0kb", "5kb", "10kb"))

pdf("Plots/KIRP_circular_correlation_meth_exp_heatmap_methreg.pdf")
plot.new()
circle.size <- unit(1, "snpc") # snpc unit gives you a square region
pushViewport(viewport(x = 0, y = 0.5, width = circle_size, height = circle_size,
                      just = c("left", "center")))
par(omi = gridOMI(), new = TRUE)
circlize_plot()

追溯：

Error in circos.link(df$rn[k], c(df$x1[k] - abs(df$value1[k]), df$x1[k]), :
unused arguments (split = c("2a", "1a", "1b", "1a", "1a", "1b", "2b", "2b", "2a", "1a", "2b", "1b", "1c", "2b", "2a", "2b", "1b", "2b", "2b", "1c", "1c", "2b", "1c", "2c", "1a", "1a", "2b", "1c", "1a", "1c", "2b", "2a", "2c", "2a", "1a", "1a", "1b", "1a", "1b", "1a", "2a", "1a", "1c", "1a", "1c", "2b", "1a", "2a", "1c", "2a", "1a", "2b", "2b", "2c", "2b", "2c", "2b", "2a", "2b", "2a", "2b", "2a", "2a", "2a", "2b", "2b", "2b", "2c", "2c", "2b", "2b", "2b", "2a", "2b", "1c", "2b", "2a", "2b", "2b", "2a", "2b",
"2b", "1a", "1b", "1c", "1a", "1a", "2a", "1c", "1a", "1a", "1c", "2a", "2b", "2a", "1c", "2a", "1a", "1a", "1c", "1a", "2c", "1a", "1b", "2c", "2a", "2c", "2b", "1a", "2c", "1a", "1a", "1c", "2b", "1a", "1c", "1b", "1c", "2c", "1c", "2b", "1a", "2c", "1c", "2a", "1c", "1c", "2a", "2b", "2a", "1c", "2c", "2a", "1c", "1b", "1c", "2a", "1a", "1a", "2b", "2c", "2a", "1a", "2b", "1c", "2a", "1a", "2c", "2c", "2a", "1b", "1b", "1b", "2b", "1a", "2b", "2b", "1a", "2c", "2c", "2

剩余代码：

upViewport()
h = dev.size()[2]
lgd.list = packLegend(lgd.subtype, lgd.meth, lgd.expr, lgd.dnam.effect, lgd.tf.role, lgd.dnam.pval, lgd.TF.pval, lgd.synergistic.dnam.tf.pval, lgd.dist,
                      max_height=unit(0.9*h, "inch"))
draw(lgd.list, x=circle.size, just="left")
dev.off()

输入：

> clin.info$subtype[1:10]
 [1] "2a" "1a" "1b" "1a" "1a" "1b" "2b" "2b" "2a" "1a"

> dput(meth[1:10,1:10])
structure(c(0.970559601230577, 0.932019744622219, 0.932425186892924, 
0.457856431008969, 0.627465074882506, 0.936638025996916, 0.783963744999747, 
0.692288921425085, 0.159672532714498, 0.916221348761207, 0.966526240990762, 
0.658133739424885, 0.934725644334937, 0.208926370054914, 0.675330760134758, 
0.947481095266047, 0.68095522917082, 0.870115202168511, 0.182194260105145, 
0.947827551572278, 0.9547414276773, 0.8870658590815, 0.930339362741965, 
0.251912133635195, 0.668048974854262, 0.929015733626052, 0.81439755126243, 
0.863820816864295, 0.262184639687776, 0.919137183610684, 0.962944267941447, 
0.816969606674613, 0.920324347738064, 0.269181280099996, 0.646018376023491, 
0.94422024213849, 0.825531134571556, 0.850541845120231, 0.383124627310674, 
0.894046817443657, 0.949954692914802, 0.471749946722121, 0.807438284971849, 
0.198210791165781, 0.419648629655657, 0.925981678207389, 0.747238112656689, 
0.479487595265121, 0.135817171373483, 0.915714899098325, 0.961506276150628, 
0.497144423630466, 0.842842269651223, 0.192980115277766, 0.34694235149816, 
0.938845955692942, 0.596872875594833, 0.829311979770309, 0.095497663882203, 
0.760281835094803, 0.958608460528636, 0.391542162408205, 0.788032332006273, 
0.0858063904427799, 0.389182161991806, 0.929826820318033, 0.597054145633082, 
0.887124132613724, 0.105961883042841, 0.889848592622701, 0.959619815668203, 
0.641986691253107, 0.852885768110978, 0.166605523819506, 0.388712417817833, 
0.909146129816036, 0.655313654779395, 0.812348993288591, 0.155677269854645, 
0.870826675637036, 0.939628877812377, 0.867474877846409, 0.901774530271399, 
0.135383203304268, 0.529803404308775, 0.908866657530492, 0.287143254078514, 
0.813718332481565, 0.235871797160968, 0.903798027381128, 0.975728903408741, 
0.449572131549706, 0.818444637178375, 0.0740835421033795, 0.344167810233288, 
0.947385012179395, 0.741077844311377, 0.763177249325815, 0.174783683559951, 
0.925959780621572), dim = c(10L, 10L), dimnames = list(c("cg03477043", 
"cg00926926", "cg00488747", "cg04452095", "cg04658243", "cg07156990", 
"cg01520454", "cg04247152", "cg02495310", "cg04578010"), c("TCGA-2K-A9WE-01", 
"TCGA-2Z-A9J1-01", "TCGA-2Z-A9JI-01", "TCGA-2Z-A9JJ-01", "TCGA-4A-A93W-01", 
"TCGA-4A-A93X-01", "TCGA-5P-A9JU-01", "TCGA-5P-A9JY-01", "TCGA-A4-7583-01", 
"TCGA-A4-7584-01")))

> dput(distal.dmr.cpgs[1:10,])
structure(list(regionID = c("chr16:2762422-2762423", "chr16:2762422-2762423", 
"chr16:89379025-89379026", "chr16:89379025-89379026", "chr16:89379025-89379026", 
"chr16:89379025-89379026", "chr16:89379095-89379096", "chr16:89379095-89379096", 
"chr16:89379095-89379096", "chr16:89379095-89379096"), probeID = c("cg03477043", 
"cg03477043", "cg00926926", "cg00926926", "cg00926926", "cg00926926", 
"cg00488747", "cg00488747", "cg00488747", "cg00488747"), target = c("ENSG00000167978", 
"ENSG00000167978", "ENSG00000167522", "ENSG00000167522", "ENSG00000167522", 
"ENSG00000167522", "ENSG00000167522", "ENSG00000167522", "ENSG00000167522", 
"ENSG00000167522"), target_symbol = c("SRRM2", "SRRM2", "ANKRD11", 
"ANKRD11", "ANKRD11", "ANKRD11", "ANKRD11", "ANKRD11", "ANKRD11", 
"ANKRD11"), TF_symbol = c("CTCF", "CREB3", "RELA", "ZNF384", 
"NFYC", "FOXK2", "RELA", "STAT3", "ZNF384", "NFYC"), TF = c("ENSG00000102974", 
"ENSG00000107175", "ENSG00000173039", "ENSG00000126746", "ENSG00000066136", 
"ENSG00000141568", "ENSG00000173039", "ENSG00000168610", "ENSG00000126746", 
"ENSG00000066136"), distance_region_target_tss = c(9795, 9795, 
111534, 111534, 111534, 111534, 111464, 111464, 111464, 111464
), target_region = c("chr16:2752626-2772538", "chr16:2752626-2772538", 
"chr16:89267630-89490561", "chr16:89267630-89490561", "chr16:89267630-89490561", 
"chr16:89267630-89490561", "chr16:89267630-89490561", "chr16:89267630-89490561", 
"chr16:89267630-89490561", "chr16:89267630-89490561"), met.IQR = c(0, 
0, 0, 0, 0, 0, 0, 0, 0, 0), RLM_DNAmGroup_pvalue = c(0.0422464065351595, 
0.0344133919449146, 0.0473017788132601, 0.0484196210021861, 0.0435840226150332, 
0.0308252756419765, 0.00150201162815966, 0.00151985887265838, 
0.00275663535734938, 0.0012627667435956), RLM_DNAmGroup_fdr = c(0.0422464065351595, 
0.0422464065351595, 0.0484196210021861, 0.0484196210021861, 0.0484196210021861, 
0.0484196210021861, 0.00189982359082297, 0.00189982359082297, 
0.00275663535734938, 0.00189982359082297), RLM_TF_pvalue = c(0.0419750640980108, 
0.0103671579957865, 0.467132225063278, 0.72093859766613, 0.272092275972533, 
0.000514794320566114, 0.891428521530259, 0.853872726027603, 0.161356975335507, 
0.823791465863571), RLM_TF_fdr = c(0.0419750640980108, 0.0207343159915729, 
0.622842966751037, 0.72093859766613, 0.544184551945066, 0.00205917728226446, 
0.891428521530259, 0.891428521530259, 0.403392438338768, 0.891428521530259
), `RLM_DNAmGroup:TF_pvalue` = c(0.816616675064483, 0.315518615926782, 
0.334800389930188, 0.886558336903593, 0.288111456068547, 0.0553116797875954, 
0.674629382494593, 0.611101293078096, 0.82494221468516, 0.601335599655604
), `RLM_DNAmGroup:TF_fdr` = c(0.816616675064483, 0.631037231853565, 
0.446400519906918, 0.886558336903593, 0.446400519906918, 0.221246719150382, 
0.82494221468516, 0.82494221468516, 0.82494221468516, 0.82494221468516
), RLM_DNAmGroup_estimate = c(0.391366841432446, 0.371327086133833, 
0.395884100676198, 0.390536544990631, 0.40695391264101, 0.408579807600067, 
0.70307776267823, 0.706552565374348, 0.649042588677745, 0.721190519859083
), RLM_TF_estimate = c(0.281227589928979, -0.323039230805149, 
-0.110519428773639, 0.0548209506490028, 0.162772915125468, 0.463202502783608, 
-0.0212827073929813, -0.0295766095061422, 0.195192883961314, 
0.0369136296170289), `RLM_DNAmGroup:TF_estimate` = c(0.0476365755291261, 
-0.165727596247364, 0.203133492532222, 0.0292733545380043, -0.216456686525197, 
-0.356607795519526, 0.0847642614797746, 0.109600973249801, -0.0437870163085532, 
0.115542672949099), Model.quantile = c("Robust Linear Model", 
"Robust Linear Model", "Robust Linear Model", "Robust Linear Model", 
"Robust Linear Model", "Robust Linear Model", "Robust Linear Model", 
"Robust Linear Model", "Robust Linear Model", "Robust Linear Model"
), Target_gene_DNAm_high_vs_Target_gene_DNAm_low_wilcoxon_pvalue = c(0.0111947515949219, 
0.0111947515949219, 0.02969053900481, 0.02969053900481, 0.02969053900481, 
0.02969053900481, 0.00113717643198146, 0.00113717643198146, 0.00113717643198146, 
0.00113717643198146), TF_DNAm_high_vs_TF_DNAm_low_wilcoxon_pvalue = c(0.417943629363925, 
0.135674453768463, 0.145235627143256, 0.680264511884829, 0.113088638511403, 
0.926405872904368, 0.417943629363925, 0.560144252514798, 0.495175592310321, 
0.179309287244037), `% of target genes not expressed in DNAm_low and DNAm_high` = c("0 %", 
"0 %", "0 %", "0 %", "0 %", "0 %", "0 %", "0 %", "0 %", "0 %"
), DNAm_low_RLM_target_vs_TF_pvalue = c(0.0259364958125814, 0.00753553623999004, 
0.460325106599497, 0.692693342638731, 0.315289731003607, 0.000756745896691946, 
0.892812259055557, 0.850393817987769, 0.137294538551965, 0.813640144669246
), DNAm_low_RLM_target_vs_TF_estimate = c(0.27906988386453, -0.323013805407332, 
-0.118743186892688, 0.0640624351376578, 0.158695080628396, 0.460733657087822, 
-0.020838062194407, -0.0301913231291102, 0.206016183618269, 0.038783711513423
), DNAm_high_RLM_target_vs_TF_pvalue = c(0.0527939209600978, 
0.000106204477959615, 0.493175812986649, 0.514588360117673, 0.728135657907864, 
0.422326547532424, 0.632621266523826, 0.585265337302442, 0.302558584728815, 
0.297918126183749), DNAm_high_RLM_target_vs_TF_estimate = c(0.333610936801836, 
-0.488728802926698, 0.0951631917249828, 0.0845092162375904, -0.0460710731343046, 
0.106818566360051, 0.0621604325281235, 0.0796459855690455, 0.149478480439067, 
0.152745219994597), DNAm.effect = c("ns", "ns", "ns", "ns", "ns", 
"ns", "ns", "ns", "ns", "ns"), TF.role = c("ns", "Repressor", 
"ns", "ns", "ns", "Activator", "ns", "ns", "ns", "ns"), cpg_type = c("distal", 
"distal", "distal", "distal", "distal", "distal", "distal", "distal", 
"distal", "distal")), row.names = c(785L, 786L, 866L, 867L, 868L, 
869L, 870L, 871L, 872L, 873L), class = "data.frame")

> dput(exp[1:10,1:10])
structure(c(8.8764930371619, 7.10418439636105, 5.82248600796462, 
9.6336088827523, 7.12980348328146, 5.93104377991884, 7.21600791720845, 
11.0446571390993, 9.43701543521221, 6.93486223041719, 8.50488041264077, 
7.4053080406626, 6.10447768176044, 9.82618201799457, 7.25526533976169, 
5.62784390424282, 7.55732162812635, 10.8857274881571, 9.21597835642978, 
7.32131574109185, 8.45739499869935, 7.02069419722442, 5.77970193951746, 
9.43056453489455, 7.16074467184323, 5.54875953204378, 7.46818856378399, 
10.2099696597136, 9.89181133862426, 6.79730741996895, 9.11863785771178, 
7.63223080712203, 5.88771420376909, 9.28640690000328, 7.34612983819779, 
6.08474958802744, 7.47638114657282, 10.5659651734064, 9.64444452544193, 
6.87091467703423, 9.73111037922987, 7.4250931669883, 6.17996170442533, 
10.1043913365228, 6.80393840392927, 6.04973861895014, 7.03770927524486, 
9.6918746764765, 9.21727922543568, 6.7071916389237, 8.31760092734971, 
7.37678634610584, 5.77890345814874, 10.2596698559209, 7.07928303958814, 
5.46545945570634, 6.26016743442868, 9.84488278708719, 9.43848400972158, 
7.12912726034445, 8.05465598032384, 7.1768739722416, 5.98202385467534, 
9.02483580988576, 7.49478518039372, 5.29219150981698, 7.53579328045101, 
11.0222325979263, 10.6179427577002, 7.84706691501966, 8.79411208840637, 
7.26461905608631, 6.3968834641627, 8.8518016637199, 7.30708053720346, 
5.08512177701773, 6.83888005212149, 10.8453676226894, 9.79866023188757, 
7.43806001401501, 8.7026874548283, 7.19605801821642, 5.84149558008759, 
10.0241110621631, 7.04542218287691, 5.68002581559096, 7.10416270371368, 
9.77279638376636, 9.31061956515004, 7.1075499703718, 8.31733643440423, 
7.20641252920776, 6.32798187568294, 9.25596850084915, 6.92549753101381, 
5.47403519965473, 6.90283808416669, 10.6081725987497, 8.94455634675817, 
7.08053348786898), dim = c(10L, 10L), dimnames = list(c("SRRM2", 
"ANKRD11", "RPTOR", "HSP90AA1", "RER1", "ZC3H18", "PDAP1", "CD63", 
"RPL5", "MAPKAPK2"), c("TCGA-2K-A9WE-01", "TCGA-2Z-A9J1-01", 
"TCGA-4A-A93W-01", "TCGA-5P-A9JY-01", "TCGA-A4-7584-01", "TCGA-A4-7585-01", 
"TCGA-AL-7173-01", "TCGA-B1-A47O-01", "TCGA-B1-A657-01", "TCGA-B9-A5W8-01"
)))

Answer 1

我认为您发布了错误的代码，或者您正在尝试使用

chordDiagram()

创建热图。

您说您想创建一个圆形热图，但我看到您发布了和弦图的代码。对于圆形热图，您应该使用

circos.heatmap()

。您可以在circlize_book

上查看编码

关于错误，则生成该错误是因为

split

参数不可用于

chordDiagram()

。您不能在

chordDiagram()

中使用它。