meme CentriMo富集 motif匹配序列的获取

背景

当我们使用 CentriMo 进行motif 富集后,如采用下面的命令

$ centrimo -oc test test.500l.fa  $meme_motif

会在test目录产生三个文件(centrimo.htmlcentrimo.tsvsite_counts.txt),用浏览器打开centrimo.html后,点击motif id前面的复选框后,会在左边位置的文本框显现出motif匹配的序列。
meme CentriMo富集 motif匹配序列的获取

因此,问题是我想要获取每个motif匹配的序列。一开始,我想应该centrimo.tsv,会有这吧,但是没找到。然后又去centrimo 的参数说明里看看有没相关参数设置,找到了个,不把匹配序列放进html的参数。。。
meme CentriMo富集 motif匹配序列的获取

既然如此那就只好自己弄了。

方法

既然能在html上显示,那网页源码里肯定包含了这些信息。所以查看下网页源码:
meme CentriMo富集 motif匹配序列的获取
分析下, 通过var data定义个JavaScript对象(相当于Python中的字典),sequences储存着所有的序列,motifs下存储着motif信息,每个motif以js对象表示,其中每个motif对象里,seqs又存储着序列索引(sequences数组的索引,没在图片中给出,可自行查看自己的的网页源码)。

因此,我们只要提取出 data 对象中的sequencesmotifs数据, 根据motifs下的每个motif储存的 seqs索引在sequences中取出对应的序列就可以获得了。

所以上代码(Python):

import re
import json
from bs4 import BeautifulSoup  ## 安装 pip install beautifulsoup4


file = "centrimo.html"

## 读取html文件
with open(file, "r") as f:
	html = f.read()


## 通过bs4 进行解析
soup = BeautifulSoup(html)

# 在script标签下,找到包含 var data的文本, \s*是用于匹配空格的
script = soup.find('script', text=re.compile('var\s*data\s*=\s*'))

## 通过正则匹配到 var data = {...}; 中data对象里的数据 
json_text = re.search('var\s*data\s*=\s*({.*?})\s*;$', script.string, flags=re.DOTALL | re.MULTILINE).group(1)

data = json.loads(json_text)  #加载成字典 
print(data.keys())

运行后,输出是这样,代表data数据提取成功

dict_keys(['version', 'revision', 'release', 'program', 'cmd', 'options', 'seqlen', 'tested', 'alphabet', 'background', 'sequence_db', 'motif_dbs', 'sequences', 'motifs'])

简单的代码就提取了。

sequences = data["sequences"]
motifs = data["motifs"]

extracted_seq = []  # 这里作为一个例子,把它们保存在列表里。

for motif in motifs:
	seqs_idx = motif["seqs"]
	seqs_id = [sequences[i] for i in seqs_idx]
	extracted_seq.append({"id": motif["id"], "seq": seqs_id})

print(extracted_seq[1]["seq"][:10])

# ['chr2:22587565-22588065', 'chrX:143482808-143483308', 'chr1:7397843-7398343', 'chr5:88764784-88765284', 'chr9:20651695-20652195', 'chr8:25016942-25017442', 'chr6:21949731-21950231', 'chr10:30842681-30843181', 'chr3:138442803-138443303', 'chr2:3283792-3284292']

之后做什么处理,就看各自的需求了。

在此,结束了。

参考

https://*.com/questions/13323976/how-to-extract-a-json-object-that-was-defined-in-a-html-page-javascript-block-us
https://meme-suite.org/meme/doc/centrimo.html?man_type=web

上一篇:Java设计模式--观察者模式到监听器


下一篇:javascript – EmberJS:无法观察@ each.property.another_property