之前看基因家族相关的文章绘制 sequences logo 基本会使用 meme+webologo 的组合,前一阵子在公众号 R语言中文社区 推送的文章中发现了R语言的 ggseqlogo 可以用来绘制 sequences logo,好多种颜色可供选择,结果非常漂亮,帮助文档十分详细,使用也相对简单。

简单记录使用meme鉴定保守基序,简单处理meme结果文件成R语言可读取的格式,最后用ggseqlogo绘图的过程

先上结果图

使用R package ggseqlogo 简单绘制sequences logo_r语言

本地meme所用到的的命令

/usr/yanming/Bioinformatic/MEME/meme_4.12.0/bin/meme PgLHC_candidate_16.fasta -protein -oc PgLHC_motif 

-nostatus -mod anr -nmotifs 5 -minw 6 -maxw 50 -maxsize 600000000

-oc 后接输出的结果文件

-mod 后接的应该是方法,有三个参数可选 oops zoops anr 具体是什么意思自己也不太懂,好像是anr参数比较好,但是相对运行时间会比较长

-nmotif 后接要搜索的motif的数量

-minw -maxw 后接motif的最短和最长的长度

-maxsize 如果蛋白的fasta文件中序列数量较多就需要相应增到后面接的数字

使用R package ggseqlogo 简单绘制sequences logo_拖拽_02

一直都会遇到报错,因为可以输出结果,也就没有仔细探究报错的原因


使用R package ggseqlogo 简单绘制sequences logo_分隔符_03
输出的结果就是对应的motif的图片,应该可以直接使用,如果要自己绘制的话需要用到meme.txt文件


使用R package ggseqlogo 简单绘制sequences logo_拖拽_04

这是结果文件中的保守基序之一,将其复制到一个新的文件中简单处理成R可以读入的格式,我用的是linux的cut命令

cut -d ")" -f 2 | cut -d " " -f 2 > ggseq_logo_input.txt

首先用括号作为分隔符提取第二列,用竖线将上一步的结果传递给下一步作为输入,再用空格作为分隔符提取第二列,再用大于号重定向到结果文件中

使用R package ggseqlogo 简单绘制sequences logo_拖拽_05
结果


使用R package ggseqlogo 简单绘制sequences logo_拖拽_06
R绘图的代码


使用R package ggseqlogo 简单绘制sequences logo_r语言_07
使用R package ggseqlogo 简单绘制sequences logo_分隔符_08
结果就是这样似的

meme在线用法

使用R package ggseqlogo 简单绘制sequences logo_拖拽_09

结果
使用R package ggseqlogo 简单绘制sequences logo_分隔符_10

第一步

PS:今天Xftp链接服务器莫名其妙的就没有了文件的下载权限,什么原因???搞不太明白,正常拖拽上传,却不能拖拽下载