美文网首页R for statisticsData science生信绘图
跟着Science学作图:R语言ggplot2画箭头展示变量对主

跟着Science学作图:R语言ggplot2画箭头展示变量对主

作者: 小明的数据分析笔记本 | 来源:发表于2022-04-06 17:09 被阅读0次

    论文

    https://www.science.org/doi/10.1126/science.abk0989

    image.png

    最近朋友圈好多人都在转这个论文,我也找来看了看,论文研究的内容看的还是一知半解。

    论文用到的数据代码都是公开的,我们可以学习一下其中的代码

    代码链接

    https://github.com/James-S-Santangelo/glue_pc

    今天的图文重复论文中的Figure 2B

    image.png

    这个图的图注写的是The eigenvectors for environmental variables, colored according to their contribution to PC2

    这里为什么只展示对PC2的贡献暂时还不明白。主要是论文的研究内容看不明白

    本篇推文只记录画图代码了

    还是先做主成分分析

    library(readr)
    dat01<-read_csv("phenotypic-analyses/sciencefig2A.csv")
    dim(dat01)
    colnames(dat01)
    dat02<-read_csv("phenotypic-analyses/sciencefig2A_group_info.csv")
    dim(dat02)
    colnames(dat02)
    
    library(vegan)
    enviroPCA <- rda(dat01, 
                     scale = TRUE, na.action = "na.omit")
    
    eig <- enviroPCA$CA$eig
    percent_var <- eig * 100 / sum(eig)
    PC1_varEx <- round(percent_var[1], 1)  # Percent variance explained by PC1
    PC2_varEx <- round(percent_var[2], 1)  # Percent variance explained by PC2
    
    

    计算物种贡献百分比

    论文里提供的代码里放了一个参考链接
    https://stackoverflow.com/questions/50177409/how-to-calculate-species-contribution-percentages-for-vegan-rda-cca-objects

    contrib <- round(100*scores(enviroPCA, display = "sp", scaling = 0)[,2]^2, 3)
    

    生成作图数据

    library(tidyverse)
    enviroPCA_vars  <- scores(enviroPCA, display = 'species', choices = c(1, 2), scaling = 2) %>% 
      as.data.frame() %>% 
      rownames_to_column(., var = 'var') %>% 
      mutate(contrib = contrib)
    

    准备配色

    library(wesanderson)
    pal <- wes_palette("Darjeeling1", 3, type = "continuous")
    

    作图主题的一些设置

    library(ggplot2)
    ng1 <- theme(aspect.ratio=0.7,panel.background = element_blank(),
                 panel.grid.major = element_blank(),
                 panel.grid.minor = element_blank(),
                 panel.border=element_blank(),
                 axis.line.x = element_line(color="black",size=1),
                 axis.line.y = element_line(color="black",size=1),
                 axis.ticks=element_line(size = 1, color="black"),
                 axis.ticks.length=unit(0.25, 'cm'),
                 axis.text=element_text(color="black",size=15),
                 axis.title=element_text(color="black",size=1),
                 axis.title.y=element_text(vjust=2,size=17),
                 axis.title.x=element_text(vjust=0.1,size=17),
                 axis.text.x=element_text(size=15),
                 axis.text.y=element_text(size=15),
                 strip.text.x = element_text(size = 10, colour = "black",face = "bold"),
                 strip.background = element_rect(colour="black"),
                 legend.position = "top", legend.direction="vertical",
                 legend.text=element_text(size=17), legend.key = element_rect(fill = "white"),
                 legend.title = element_text(size=17),legend.key.size = unit(1.0, "cm"))
    

    最后的作图代码

    enviroPCA_variableContrib <- ggplot() +
      geom_hline(yintercept = 0, linetype = "dotted") +
      geom_vline(xintercept = 0, linetype = "dotted") +
      geom_segment(data = enviroPCA_vars, aes(x = 0, xend = PC1, y=0, yend = PC2, color = contrib), 
                   size = 2, arrow = arrow(length = unit(0.02, "npc")), alpha = 1) +
      geom_text(data = enviroPCA_vars,
                aes(x = PC1, y = PC2, label = var,
                    hjust = "inward", vjust =  0.5 * (1 - sign(PC1))),
                color = "black", size = 3.5) + 
      xlab(sprintf("PC1 (%.1f%%)", PC1_varEx)) + ylab(sprintf("PC2 (%.1f%%)", PC2_varEx)) +
      scale_colour_gradientn(colours = rev(pal), breaks = seq(from = 5, to = 25, by = 5)) +
      # scale_x_continuous(breaks = seq(from = -1, to = 1, by = 0.25)) +
      # scale_y_continuous(breaks = seq(from = -1, to = 1, by = 0.25)) +
      ng1 + theme(legend.position = "top",
                  legend.direction="horizontal",
                  # legend.title = element_blank(),
                  legend.key.size = unit(0.5, "cm"),
                  legend.spacing.x = unit(0.1, "cm"),
                  legend.text = element_text(size=10)) +
      guides(color = guide_colourbar(barwidth = 10, barheight = 0.5))
    enviroPCA_variableContrib
    
    image.png

    本期推文的示例数据和代码可以给推文点赞 点击在看 然后留言20220406获取

    欢迎大家关注我的公众号

    小明的数据分析笔记本

    小明的数据分析笔记本 公众号 主要分享:1、R语言和python做数据分析和数据可视化的简单小例子;2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记;3、生物信息学入门学习资料及自己的学习笔记!

    相关文章

      网友评论

        本文标题:跟着Science学作图:R语言ggplot2画箭头展示变量对主

        本文链接:https://www.haomeiwen.com/subject/aborsrtx.html