論文是
Whole-genome resequencing of 445 Lactuca accessions reveals the domestication history of cultivated lettuce

image.png

這篇論文的數(shù)據(jù)是公開的，代碼也公開了一部分，那我們就可以按照他的代碼來學(xué)二代測序的數(shù)據(jù)分析啦

今天我們來學(xué)習(xí)一些論文中 Fig1c 的畫圖代碼，如下

image.png

數(shù)據(jù)對應(yīng)的是論文中的 source data figure1

圖的主要內(nèi)容是散點(diǎn)圖展示主成分分析的結(jié)果，并且將局部的區(qū)域放大展示

首先是讀入數(shù)據(jù)

df<-readxl::read_excel("NG/41588_2021_831_MOESM4_ESM.xlsx",
                       sheet = "Fig1c",
                       n_max = 441)
head(df)
tail(df)

這里直接讀入excel文件用到的是readxl包中的read_excel()函數(shù)，需要制定

• 文件路徑
• excel表格中的sheet名稱

這里的n_max參數(shù)是指定讀進(jìn)來的數(shù)據(jù)的最多行數(shù)，英文這個數(shù)據(jù)集結(jié)尾處有一些注釋內(nèi)容，我們不需要，所以需要制定這個參數(shù)，自己的數(shù)據(jù)集通常是不需要指定這個參數(shù)的

先畫一個簡單的散點(diǎn)圖

library(ggplot2)
ggplot(data=df,aes(x=PC1,y=PC2))+
  geom_point(aes(color=Species))

image.png

這個看著和論文中的有些不一樣，仔細(xì)看看應(yīng)該是 論文中對 PC1去了一個負(fù)數(shù)，而且論文中的圖也映射了點(diǎn)的形狀

對PC1取一個負(fù)數(shù)

library(dplyr)
df %>% 
  mutate(PC1.1 = - PC1) -> df

畫圖

ggplot(data=df,aes(x=PC1.1,y=PC2))+
  geom_point(aes(color=Species,shape=Species))

image.png

這里會遇到一個警告信息

Warning messages: 1: The shape palette can deal with a maximum of 6 discrete values because more than 6 becomes difficult to discriminate; you have 9. Consider specifying shapes manually if you must have them. 2: Removed 369 rows containing missing values (geom_point).

形狀超過留個就難以區(qū)分了，如果你非要用超過6個的形狀，這里需要手動指定

這里是9個形狀
手動指定

ggplot(data=df,aes(x=PC1.1,y=PC2))+
  geom_point(aes(color=Species,shape=Species),
             size=5)+
  scale_shape_manual(values=15:23)

image.png

接下來是簡單的美化了

ggplot(data=df,aes(x=PC1.1,y=PC2))+
  geom_point(aes(color=Species,shape=Species),
             size=5)+
  scale_shape_manual(values=15:23)+
  theme_bw()+
  theme(legend.position = "none",
        panel.grid = element_blank())+
  labs(x="PC1 37.01% of variance",
       y="PC2 22.09% of variance")+
  geom_rect(aes(xmin=0,xmax=0.04,ymin=-0.01,ymax=0.02),
            fill="white",alpha=0,
            color="black")+
  annotate(geom="text",x=0.01,y=0.01,label="GP1")+
  annotate(geom="text",x=-0.09,y=0.05,
           label="italic(L.saligna)",parse=T,
           color="#00b9e3")+
  annotate(geom="text",x=-0.02,y=-0.075,
           label="italic(L.georgica)",parse=T,
           color="#00ba38")+
  annotate(geom="text",x=-0.07,y=-0.14,
           label="italic(L.virosa)",parse=T,
           color="#ff82cf")

image.png

那接下來是如何放大展示呢？想到了兩種方案

• 第一個是拼圖
• 第二個是借助ggforce這個包里的facet_zoom()函數(shù)
  可以參考 https://github.com/thomasp85/ggforce/pull/202

具體如何實(shí)現(xiàn)有時間再來研究吧！

歡迎大家關(guān)注我的公眾號

小明的數(shù)據(jù)分析筆記本

小明的數(shù)據(jù)分析筆記本 公眾號 主要分享：1、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡單小例子；2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)、群體遺傳學(xué)文獻(xiàn)閱讀筆記；3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！